2900- 392 - T Report of Project MICHIGAN THE ANALYSIS OF A CERTAIN CLASS OF NONLINEAR SYSTEMS A. Y. BILAL September 1963 Navigation and Guidance Laboratory Ev444&e 4 Sc4ee aW 7eEdy THE U N I V E R S I T Y OF M I C H I G A N Ann Arbor, Michigan

Institute of Science and Technology The University of Michigan Sam NOTICES Sponsorship. The work reported herein was conducted by the Institute of Science and Technology for the U. S. Army Electronics Command under Project MICHIGAN, Contract DA-36-039 SC-78801 and for the National Science Foundation under grant GP-524. Contracts and grants to The University of Michigan for the support of sponsored research by the Institute of Science and Technology are administered through the Office of the Vice-President for Research. Note. The views expressed herein are those of Project MICHIGAN and have not been approved by the Department of the Army. Distribution. Initial distribution is indicated at the end of this document. Distribution control of Project MICHIGAN documents has been delegated by the U. S. Army Electronics Command to the office named below. Please address correspondence concerning distribution of reports to: Commanding Officer U. S. Army Liaison Group Project MICHIGAN The University of Michigan P. O. Box 618 Ann Arbor, Michigan DDC Availability. Qualified requesters may obtain copies of this document from: Defense Documentation Center Cameron Station Alexandria, Virginia Final Disposition. After this document has served its purpose, it may be destroyed. Please do not return it to the Institute of Science and Technology.

Institute of Science and Technology The University of Michigan PREFACE Project MICHIGAN is a continuing long-range research and development program for advancing the Army's combat surveillance and target-acquisition capabilities. The program is carried out by a full-time staff of specialists in physics, engineering, mathematics, and psychology at the Institute of Science and Technology, and by members of the teaching faculty and graduate students of other research groups and laboratories of The University of Michigan. The emphasis of the Project is upon research in imaging radar, MTI radar, infrared, radio location, image processing, and special investigations. Particular attention is given to all-weather, long-range, high-resolution sensory and location techniques. Project MICHIGAN was established by the U. S. Army Signal Corps at The University of Michigan in 1953 and has received continuing support from the U. S. Army. The Project constitutes a major portion of the diversified program of research conducted by the Institute of Science and Technology in order to make available to government and industry the resources of The University of Michigan and to broaden the educational opportunities for students in the scientific and engineering disciplines. Progress and results described in reports are continually reassessed by Project MICHIGAN. Comments and suggestions from readers are invited. Robert L. Hess Director Project MICHIGAN

Institute of Science and Technology The University of Michigan ACKNOWLEDGMENTS The assistance and encouragement received from the administration of The University of Michigan, the faculty of the Department of Electrical Engineering, and the staff of the Navigation and Guidance Laboratory of the University's Institute of Science and Technology were instrumental in the completion of this work. The interest and continual aid given by Professor L. F. Kazda and Professor W. Porter are especially acknowledged. The suggestions, comments, and criticisms of Professors L. F. Kazda, W. Brown, Kuei Chaung, W. Kaplan, and Arch Naylor are sincerely appreciated. The author also wishes to express his thanks to Mrs. Louise Snider, who edited this report, Mr. August Antones, who drew the figures, and Mrs. Linda Bennett, who typed the first and final drafts.

Institute of Science and Technology The University of Michigan CONTENTS Notices............................. ii Preface.................................... iii Acknowledgments................................. iv List of Figures................................ vii List of Tables...........o.....................viii List of Symbols................................ ix Abstract..................................... 1 1. Introduction................................. 1 1.1. Historical Origin and Review 1 1.2. Summary and Research Objectives 20 2. A Theory for the Analysis of a Certain Class of Nonlinear TimeVarying Physical Systems........................ 20 2.1. System Description and Definitions 21 2.2. Mathematical Restrictions and the System Analysis 22 2.3. Theorem A 26 2.4. Theorem B 29 2.5. Theorem C 30 2.6. The Approximation of the Solution in the Mean Square Sense 32 2.7. Calculation of an Upper-Bound Function P(t) for the System Response x(t) 36 2.8. Variation of the Solution with Respect to the Initial Conditions and the System Parameters (Theorem D) 37 2.9. Three Example Calculations 40 2.10. Conclusions and Remarks 52 3. Stability Analysis........................ 55 3.1. System Restrictions 55 3.2. The State-Variable Representation 56 3.3. The Properties of the State Variables (v1, v2,..., vn) of the Class of Systems Under Consideration (Theorem E) 58 3.4. Stability Considerations 62 3.5. An Example Calculation: The Analysis of a Second-Order Nonlinear Control System 67 3.6. Conclusions and Remarks 76 4. Extension of the Analysis to Multiloop Systems 77 4.1. Introduction 77 4.2. System Description and Definitions 78 4.3. The Mathematical Restrictions and the System Analysis (Theorems F, G, and H) 83 4.4. Phase Space Representation and Stability Analysis (Theorem I) 91 4.5. Example Calculations 95 4.6. Conclusions and Remarks 101 5. System Analysis by Numerical Techniques 102 5.1. Introduction 103 5.2. The Estimation of Trucation Errors by Two Measurements 104 5.3. The Systematic Approximation to the Partition Method 106

Institute of Science and Technology The University of Michigan 5.4. The Accuracy of the Approximate Solution 109 5.5. An Example Calculation Using the Digital Computer 115 5.6. Comparison with Similar Methods 122 5.7. Conclusions and Remarks 123 6. Conclusions, Criticisms, and Suggestions for Future Study......... 124 6.1. Summary 124 6.2. Questions for Future Work 126 Appendix A. The Exact Solution of Equation 30................ 128 Appendix B. A Summary of the L2 Space.................... 129 Appendix C. The Impulse Response of Time-Varying Linear Systems and the Extension of the Analysis of Section 2 to Situations in Which the Forcing Function Is gl(t) = K(D, t)g(t) Where K(D, t) Is a Linear Time-Varying Operator......... 130 Appendix D. A Proof for Equation 133...................... 132 Appendix E. Uniqueness and Successive Approximations........... 133 Appendix F. Extension of the Analysis to More Complex Nonlinearities........................ 134 Appendix G. A Proof of Equation 289......... 135 Appendix H. A Continuous Function Which Belongs to L2[0, co] and Does Not Have a Limit at t = +oo................... 137 Appendix I. Convergence in the Mean and Convergence in the Ordinary Sense....................... 138 Appendix J. Interpolation, Approximate Differentiation, and Approximate Integration.................... 139 References................................... 149 Distribution List...................... 152 Vl

Institute of Science and Technology The University of Michigan FIGURES 1. Typical Outline of Steps Followed in Finding an Approximate Solution for a Nonlinear Equation... 3 2. The Physical Meaning of x0(t)......................... 6 3. The Physical Meanings of x1(t) and M0(s)................... 6 4. The Physical Meanings of x2(t) and M1(s)................... 6 5. The Physical Meanings of Xn(t) and M 1(s)................... 7 6. The Partitioned Nonlinear System....................... 8 2n 7. A Plot of Versus n........................... 36 (n- 1V s 8. An RL Circuit with a Nonlinear Resistor................... 40 9. Equivalent Representations of the System Governed by Equation 185.....41 10. A Plot of the Solution x(t) and the Upper-Bound Function P(t) for the System Governed by Equation 196.....................45 11. Representation of the System Governed by Equation 221........... 45 12. A Plot of d(t) and:(t) for the System Governed by Equation 237........ 50 13. Electric Circuit Configuration of Equation 96...................54 14. A Feedback System Representation of the Nonlinear Integral Equation 122..54 15. A Linear System Representation of Equation 97.............. 54 16. System Trajectory............................... 62 17. Projection of a Solution Curve......................63 18. Geometric Definition of Stability........................64 19. A Second-Order Nonlinear Control System..................68 20. Required System Characteristics....................... 70 21. The System Transient Response....................75 22. Phase Plane Trajectory of the Second-Order Nonlinear Control System..76 23. Actual and Average System Behavior......................77 24. A Two-Terminal Pair, Electric Circuit Representation of the Differential Equation 345............................80 25. A Multiloop Nonlinear Feedback System.................... 81 26. The Conjugate Multiloop Nonlinear Feedback System.............81 27. A multiloop Nonlinear Feedback System with Four Nonlinear Elements, Two Linear Elements, and Five Amplifier Gains............... 96 28. A Multiloop Nonlinear Circuit Equivalent to that of Figure 27......... 97 29. A Two-Terminal Pair Network with Nonlinear Resistors and Nonlinear Coupling....................................97 30. The Measurement of x1......................... 106 vii

Institute of Science and Technology The University of Michigan 31. A Plot of 6 and $ Versus h.................... 111 r n 32. A Plot of the Actual Error Propagation and the Theoretical Error Propagation for Different Values of A 0.119 33. A Typical Flow Chart for the Solution by the Method of Successive Approximation............................... 120 34. Different Approximation for the Actual Solution.............. 122 35. A Continuous Function f(t) Which Belongs to L2[10, c] and Does Not Have a Limit at t = +....................... 137 36. Geometrical Meaning of the Weierstrass Theorem............. 139 TABLES I. Some Numerical Values of x(t) i and P(t) for the System Governed by Equation 196..............................45 II. Some Numerical Values of d(t) and P(t) for the System Described in Equation 237................51 III. Some Numerical Values of vl(t), v2(t) for the Second-Order Nonlinear Control System................... 75 IV. The Mean Square Error q2 for Different Values of Iterations n.......116 V. Numerical Results from the Digital Computer...............121 VI. Diagonal Difference Table..................140 VII. Horizontal Difference Table....................... 141 VIII. Central Difference Table...........................143 * * Vlll

Institute of Science and Technology The University of Michigan SYMBOLS ai (t) = System time-varying parameters 2 = A number obtained from T 2 2 A = A number obtained from a (t) dt < A < 0 T 2 C = A number obtained from n 2(t)dt < m - m 0 C2 = A number obtained fromT n (t)dt < Cm2 n n Dn, the n-th differential operator n = 0, 1, 2,... p dtn d(t) = An enforced bound within which the system response must remain E(i k) = ET[X (i)(kh)], the truncation error in the evaluation of xm (kh) m E = Round-off errors r ET(xt) = The truncation error in x(t) F(t, u, x) = The weighted nonlinear function (for the single-degree-of-freedom case) = W(t, u)N(x) Fik(t, u, x) = Wii(t, u)Nik(x) = the weighted nonlinear functions (for multiloop case) F(m)(t, u, x) = The m-th partial derivative of F(t, u, x) with respect to t F. (t, u, z)= The /i-th derivative of Fik(t, u, z) with respect to time ik ik gl'.. g = System inputs for the multiloop system g(t) = System input for the single-degree-of-freedom case h = Subinterval of time total interval number of intervals K(t, u) = An upper bound of aF(t, u, x) in a given domain D ax Kik (t, u)= An upper bound of in a given domain D ~~~ik a~~z K (t, u) = An upper bound of aF (m)(t u x)in a given domain D KMt, ax inagiendomainax L2[0, T] = L2 space in the closed interval [0, T] L(D, t) = Time-varying linear differential operator (single-degree-of-freedom case) Lik(D, t) = Time-varying linear differential operator (multiloop case) Lim = Limit in the ordinary sense ix

Institute of Science and Technology The University of Michigan L.i.m. = Limit in the mean square sense n(t) = A function obtained from F(t, u, x0) dt < n(t) 0 nik(t) = A function of time obtained from 5 Fik()(t, u, xk0) du < nik(t) 0 n (t) = A function obtained from F(m)(t, u, x) du < n (t) m 00 m N(x) = Nonlinear functions with memory (single-degree-of-freedom case) Nik(x) = Nonlinear functions with memory (multiloop case) 2.T 2 q = Mean square error = (x -x ) dt n t = Time T = Total interval of time v(t) = State vector with component v1(t),..., v (t) representing the state variables of the system W(t, u) = Impulse response of L(D, t) Wii(t, u) = Impulse response of Lii(D, t) ixll = The L2 norm of x = x (t)dt 0 x1, x2,.., xn = System outputs for the multiloop system x(t) = System output (for the single-degree-of-freedom case) xi (t) = Solution of Lii(D, t)x(t) = gi(t) xk(i) = An approximate value of xm(i)(t) at a time t = kh'~mk m x (i)(t) = The i-th derivative of x(t) xm(t) = The m-th iterate of x(t) {xn(t)} = An infinite sequence of functions x (t), n = 1, 2,..., oo xo(t) = Solution of L(D, t)x(t) = g(t) 2 2 a (t) = A function obtained from K (t, u)du < a (t) 0 2 C 2 (t a (t) = A function obtained from 0 Km(t, u) du < a (t) 2(t) = A function of time obtained from Kik (t,u) du (t) f(t) = An upper bound for x(t) inside d(t) which belongs to L2[0, T] i(t) = Upper bounds for the state variables vi(t). ji(t) belongs to L2[0, T] x~~~~~~~~~~~~~~~~~~

Institute of Science and Technology The Univer'sity of Michigan 6(t - u) = The unit impulse response applied at t = u e = Belongs to en = The error in x(t) when approximated by x (t) = x - x (t) n n n W' (i) (i) "(i, k) = The discretization error in xmk = xmk - x (i)(kh) (I~k) mk mk mXl

THE ANALYSIS OF A CERTAIN CLASS OF NONLINEAR SYSTEMS ABSTRACT In nonlinear analysis, the partitioning technique has been used to analyze a certain class of nonlinear systems whose dynamic behavior can be represented by a nonlinear differential equation with time-varying parameters. When suitable restrictions were placed on the linear, nonlinear, and forcing function terms, the system equation presented a unique solution which existed to the right of the initial state. The system solution was given as a limit of a sequence of Picard iterates {xn} which are well defined in a given domain and which belong to L2 space. A formula was developed which permits determining the number of iterates necessary for the approximation of the solution in the mean square sense. Within the restrictions imposed on the system, the system response was found to be uniformly continuous with respect to the initial conditions and the system parameters. Two different definitions of the norms were selected, and the behavior of the system trajectory was investigated for two important cases: (1) at any instant of time and (2) on the average during an interval of interest. The system (under the given conditions) was asymptotically stable in the sense of Lyapunov. The analysis was extended to include systems whose behavior is governed by simultaneous nonlinear differential equations. It was found that a class of multiloop nonlinear systems exists for which it is possible to find the exact analytic solution to any degree of accuracy. For situations in which the forcing function, the nonlinear functions, and the time-varying parameters can be defined graphically or tabularly, a systematic method has been presented to include numerical analysis. The method lends itself easily to digital computer calculations. INTRODUCTION In literature dealing with the general area of automatic control, information on the use of the partitioning technique for the analysis and synthesis of a class of nonlinear systems has appeared only recently. This introduction is intended to provide adequate background information on this subject for the reader. The origin of the partitioning technique is discussed from the author's point of view, and this is followed by a brief review of the theory and its limitations as presented by A. A. Wolf [41 through 47]. The author also considers questions and arguments that have been raised recently as a result of comparison of the partition method with other methods of analysis. Finally, the objectives of this study are summarized and discussed. 1.1. HISTORICAL ORIGIN AND REVIEW During the past two decades, many papers have been written about techniques for analyzing and synthesizing nonlinear systems. But the techniques thus far developed are restrictive be

Institute of Science and Technology The University of Michigan cause each applies to a particular class of problems, each has its own particular value and usefulness, and each lacks an underlying or unifying theory that would join it to others. However, there are some exceptions. One is the general theory of stability via Lyapunov's Second Method [14]. This has recently attracted the attention of the engineering community of this country. Another exception is the technique of optimum control analysis. Generally, some of the techniques for analysis of nonlinear systems are linearization techniques and others are approximation techniques, sometimes quite restricted. However, the purpose of this report is not to discuss different possible approaches to the solution of nonlinear systems (since these approaches are adequately treated elsewhere [21, 44]), but to present a detailed study of one approach, the partitioning technique. This technique is applicable to a class of nonlinear systems which contain some linear terms, some nonlinear terms, and a forcing function. However, the idea of separating the original nonlinear equation into linear and nonlinear terms is not original. It was first introduced by Henri Poincare [33] at the end of the nineteenth century. In his classical researches on celestial mechanics, Poincare' developed quantitative methods of approximations by expansion in terms of small parameters. During the past thirty years or so, the methods of Poincare have been adapted to the problems of applied science in general. This was done as a result of the work of Mandelstom, Papalexi, Andronov, Krylov, Bogolyubov, and other Russian scientists working jointly in this field [15, 19, 28, 33]. In a recent book [10], Cunningham mentioned the idea of partitioning in an explanation of procedures for finding the approximate analytical solution for a nonlinear equation. The procedures are shown in Figure 1. In this diagram the original nonlinear differential equation is first separated into two parts: one part is a linear equation simple enough to permit exact solution; the other part contains any terms that are difficult to handle and usually involves the nonlinear terms, but there may be other terms as well. Next, the linear equation is solved to give the zero-order or generating solution. This generating solution is then employed in some way with the nonlinear terms of the original equation to produce first-order correction terms. Finally, these correction terms are combined with a generating solution to yield a first-order corrected solution which is an approximate solution of the original equation. The exact form of the correction terms depends upon the particular details of the method being employed. In some cases, the correction terms are merely added to the generating solution. In other cases, the correction terms produce changes in amplitude and phase of the generating solution. If the degree of nonlinearity is sufficiently small, a single application of this method will probably yield a first-order correction solution of sufficient accuracy. If the degree of nonlinearity is large, it is sometimes possible to improve the accuracy by applying the method a second time to obtain a second-order corrected solution (as shown in Figure 1). In theory, continued increases in accuracy are possible by further repeated applications of the method, but this is

Institute of Science and Technology The University of Michigan Zero-Order First-Order Linear Generating Corrected Equation o tion Solution Solution Second- Order FIGURE 1. TYPICAL OUTLINE OF STEPS FOLLOWED IN FINDING AN APPROXIMATE SOLUTION FOR A NONLINEAR EQUATION usually not practical because any small increase in accuracy is accompanied by a large increase in mathematical complications. Actually, an important uncertainty inherent in this kind of method is the error in the solution it yields, since this error is often difficult to calculate. Some well-known methods of analysis depend on the above procedure, such as the perturbation, reversion, variation of parameter, and averaging techniques by Galerkin and Ritz. Although details of the different methods vary, most of them are rather similar and follow the steps diagrammed in Figure 1. However, it should be noted that uniqueness of solution was assumed in the previous discussion. In his partition theory, A. A. Wolf studied conditions which must be satisfied by the linear, nonlinear, and forcing function terms of the system equation if the system is to yield a unique solution. Although nonlinear differential equations generally are known not to possess unique solutions [47, Appendix A, B, C, E], it is apparent from physical conditions that, for a particular initial state, physical systems must possess unique solutions which exist in a domain to the right of that initial state. This means that the solution of the given nonlinear differential equation corresponding to a given initial state should be unique; otherwise the analysis given does not adequately describe the behavior of the system Since the publication of Wolf's dissertation [47], which gave a detailed account of the partition theory, many articles have appeared which expound the deterministic process of analysis soltios wichexit i a oman t th riht f hatiniialstae. hismeas tat he olui 3

Institute of Science and Technology The University of Michigan of a class of nonlinear systems. With the aid of the statistical transform theorem, the theory of systems subjected to stochastic processes is extended [43, 47]. Let us now consider a typical class of physical nonlinear systems whose behavior can be described by the following nonlinear equation:1 Z(D)x(t) + N(x,,..., x())= g(t) (1) where the following assumptions and restrictions, which are sufficient but not necessary for the analysis, are made: (1) The operator Z(D) is linear with constant coefficients and is restricted so that it has an impulse response y(t) defined as a1 +Joo st y(t) = j J Z(s)ds (2) which is of exponential type and order at most2 one. (2) The forcing function g(t) is of exponential type and order at most one. (3) (3) The nonlinear function N(s) is single valued, continuous, and satisfies the following condition: Given a pair of positive constants M and a and two continuous functions v(t) and u(t), which -at are asymptotically like e, then N must satisfy the modified Lipschitz condition: IN[v(t)]- N[u(t)] <_ Me at Iv(t) - u(t) (4) It should be noted that condition 4 is true only when N is a function of x alone, and does not depend on x, x, or x(m). As a result of imposing these three conditions on the system described in Equation 1 and applying the partitioning techniques described later in this section, a system solution x(t) is found which is analytic and unique. The exact solution is represented as a limit of a sequence 1The method of partitioning discussed in this section is applicable only to systems described by equations of the same type as Equation 1. A function f(t) is of exponential type if it satisfies the inequality If(t) < Me at, where M and a are two arbitrary finite numbers. See Reference 37, page 248, for a precise definition of the order of an analytic function. For a discussion of situations in which N contains derivatives of x, see Reference 47, page 22. 4~~~~~~~~~~~ -1 _ 0r 3-C'41_

Institute of Science and Technology The University of Michigan of iterates {xn} as n tends to infinity. Each of these iterates {xn} (for all finite n) is of exponential type and order one, and is given by the following recurrent equation: t Xn+1 x0(t)- J N[Xn(T)] y(t - T)dT, (n = 0, 1, 2,...) (5) where, for the case of zero initial condition, x0(t) is defined by xo(t) =i g(T) y(t - T) dT (6) Therefore, x(t) = Lim x (t), n > 1 (7) n-co n This limit has been shown [47] to be a uniquely convergent procedure which gives the required solution. In other words, there is a certain class of linear systems (not necessarily physically realizable) in which the response x (t), for n = 0, 1, 2,... of each of these systems forms a point set, where the limit as n - oo is the response of the nonlinear differential Equation 1. The members of the set are generated according to Equations 5 and 6. Thus, among the sequence {Xn}, there is a particular x whose behavior approximates that of the actual response as closely -xn n as possible. This means that the iterates xn for all n > q approximate the solution x(t) in such a way that the error en as given by Ix(t)- Xn(t) = en(t) (8) satisfies the inequality E_(t) <6 (9) where 6 and q are positive numbers. The physical meaning of the generation of the iterates is explained in Figures 2 through 5. From these figures one can consider x0(t) to be the response of a linear network having zero feedback, as shown in Figure 2. One may also consider the sequence of functions x as being generated by a collection of feedback systems, the response of one member of which depends on the response of the one just preceding it. The networks Mn(s) have the property that if xn (t) excites M1 (s), the response is N[xn l(t)]. In this sense one may consider the networks Mn(s) as linear, although they may not be physically realizable. It should be noted that the method of analysis introduced above is similar in spirit to the techniques used in the so-called fixed point theorems in functional analysis. Detailed discussions of this appear in References 3 and 16.

Institute of Science and Technology The University of Michigan g(t) Y(s) x(t) FIGURE 2. THE PHYSICAL MEANING OF x0(t) g(t) g (s)x) W0(s) z }(S (a) x1(t) W M0(S) N[x0(t)] (b) FIGURE 3. THE PHYSICAL MEANINGS OF xl(t) AND M0(s) (a) xl(t). (b) M0(s). g(t) Y(s ) = (a) x2(t) M1(S) N[1( )] (b) FIGURE 4. THE PHYSICAL MEANINGS OF x2(t) AND Mi(s). (a) x2(t). (b) Ml(S).

Institute of Science and Technology The University of Michigan + 1 g(t) Y(s) = - x(t) Mn- (s) (a) x(t)' Mn 1(s) N[Snl(t)] (b) FIGURE 5. THE PHYSICAL MEANINGS OF x (t) AND M (s). (a) x (t). (b) M 1(s). The actual desired solution of the system equation can be obtained in the series form, to any degree of accuracy, by applying the rules of partitioning as follows: (1) The dynamic equation representing the system is partitioned so that the linear terms involving the highest order derivative are contained in one member of the equation. (2) A modified forcing function A(t) is then obtained which consists of the actual forcing term, the nonlinear term, and possibly some linear terms. (3) Next, the expansion of A(t) is determined by the particular properties of the linear terms, nonlinear terms, and driving function of the system, and the required form of the system solution. (This may be clarified by a consideration of the different forms of expansion of A(t) given in Equations 48 through 50.) (4) The properly partitioned differential equation is then solved. (5) Certain coefficients which arise from the auxiliary equation can then be calculated. (6) Finally, these coefficients, which appear in the partitioned equations, are eliminated to obtain the exact solution. Applying the partitioning rules to Equation 1 gives Z(D) x(t) = A(t) (10) which is the partition equation, and A(t) = g(t) - N[x(t)] (11)

Institute of Science and Technology The University of Michigan which is the auxiliary equation. The partitioned nonlinear system represents Equations 10 4 and 11 and is shown diagrammatically in Figure 6. To effectively use Equations 10 and 11 we must take advantage of the a priori properties of A(t) which were deduced by means of the conditions required of the linear and nonlinear parts of Equation 1 and its forcing function. A(t) Y(s) x(t) Linear a be O N x(t) g(t) Nonlinear N {x(t)x FIGURE 6. THE PARTITIONED NONLINEAR SYSTEM From these conditions, x(t) is deduced to be analytic. If g(t) is analytic and N is analytic, then A(t) is analytic, since it is the difference between two analytic functions. By the above means, A(t) can be expanded into a power series. In these conditions the solution of Equation 1 is given by the following: A(t)= Z C tn (12) n=O Assuming zero initial conditions for the purpose of simplicity, x(t) is given by t x(t) = f y(t - T) A(T) dT (13) t 0O t x(t)= y(t- ) C T dr=dT C Tn y(t T) )d 0 =n=0 Therefore, 0o X(t) = Z Cn Qn (14) n=O 4The configuration shown in Figure 6 is sometimes known as the canonical form of nonlinear systems.

Institute of Science and Technology The University of Michigan where Qn = the n-th moment of the impulsive response of the linear part Z(D) over the half open interval (0, t), that is, t Qn(t) = X T y(t - T)dT (15) It should be noted that if all but the highest order derivative is transposed to the right member, a power series solution results since the resulting moment function is a power series. This can be easily shown by the following reasoning: Partitioning Equation 1, as shown below, gives dxm m g(t) - M(x) (16) dtm where M(x) = N(x) + linear terms (17) for this case, therefore, we have m-l1 y(t) = (m (18) Applying Equation 15 gives Therefore, Q n(t) (m + n) t.) Equations 14 and 19 show the validity of the previous statements. Under the conditions placed on the system Equation 1 and under a special type of partitioning (at the highest-order derivative), the solution obtained is in a power series form; this fact gives rise to many questions concerning the analogy between the partition method and the classical method of analysis by Frobenius [13]. In a discussion of the partition theory [42], Professor R. L. Cosgriff (Ohio State University) stated that the partition theory could be avoided by extending the method of Frobenius to nonlinear analysis. He also stated that the present method only introduces a new facet concerning the solution of nonlinear equations, and that as each facet is exploited new insights and techniques develop, which expand our knowledge of nonlinear systems.

Institute of Science and Technology The University of Michigan Also, E. V. Bohn [5], in his comments about the partition method and associated transform methods, stated that it is difficult to agree with A. A. Wolf on the usefulness of his methods for theoretical and systematic studies, and that a power series representation is generally of little use. Bohn also pointed out that conventional methods and the partitioning method are equivalent in the following manner. The nonlinear differential Equation 1 can be put in the form un = f (u0, u1, *. * Un-1 t) (20) where the notation un (t) - d x(t) (21) dtn is used. The solution of Equation 20 can be represented by a power series co a x(t) = tn (22) n=0 From Equations 21 and 22 it follows that ak = uk(O) (23) The values of ao, al,..., an_1 can be specified as initial conditions and ak can be found by differentiating Equation 20 successively with respect to t and then substituting t = 0, which yields af af af n n n un+ u +... + u +- =f u+ 1 u n t) (24) Un+1 au0''' au n at n+1'''''n-l'Un' 0 n-i n-1 and for t = 0 an =fn (ao, al,... a n-1 O) (25) an+l = fn+l (a' al' an, ) Equation 25 should then be solved for each ak. To see the equivalence, assume that Equation 1 is partitioned at the highest derivative, as already shown by Equation 16, and let m cc -dx = E Ck t (26) 10

ititute of Science and Technology The University of Michigan From Equations 26 and 22 the following relation is obtained: aa+k Ck= ka (k>O) (27) However, the author feels that the criticisms of Bohn and Cosgriff are not totally warranted. A review of Wolf's published work would reveal that the method of Frobenius is a special case of the partition theory. The Frobenius method only gives rise to a power series solution, but the partition method gives rise to a general class of solution forms.'The nature of these solutions depends on the point of partitioning and the form of expansion of the auxiliary forcing function which, together, provide flexible control of the nature of solution forms. Reference 43 contains examples which illustrate this. They show that although the solution exists and is unique, the form of the solution is not unique but is a function of the point of partition. Thus, by using the partition technique, the solution can be obtained as a power series, Dirichlet series, trigonometric series, or orthogonal series, whereas the Frobenius method permits only a power series solution. Therefore, as a result of the partition theory, the following expansion of A(t), for example, can be justified: cc -X t A(t) = Cn e n (28) n=O and if X = an (29) n the system Equation 1 now has the following solution: 00 x(t)= Z C P (t) (30) n=O where Pn(t) are the exponential moments of the folded impulsive response of the linear part given by -nat It is evident from expanding e into a power series that a relation exists between the exponential and the impulsive moments of y(t - T), Qn(t). By choosing A as in Equation 29 we have the special case of the Dirichlet series known as the power Dirichlet series. This terminology is 11

Institute of Science and Technology The University of Michigan -at adopted since the solution is given in a power series of Z, when z = e and X is given by n Equation 29. As an example, consider a system described by [47]: dx 2 dt+ 6x+ ax =(32) where 6, o, i are constants. The initial conditions are given as x(O) = 0. Equation 32 may be partitioned in two ways: one way would be at the highest derivative and the other way would include the x itself in addition to the derivative. Consider using the highest derivative first: dx 2 (33) Thus dx C e -nat (34) dt n n=O The auxiliary function now becomes C eant x - (x2 (35) 71C e n n=O Integrating once, Equation 34 becomes Cn -nat -na e + (36) -na n=l where k = constant of integration. Substituting Equation 36 into Equation 35 results in the following: 0 00o C 0o C 0LCenat C -nat k n -nat (37) n=O n=1 n=O Let -at z e (38) so that n=O n= n=1 +2 kE n-a -k (39) n=a 12

Institute of Science and Technology The University of Michigan Upon integrating, we notice that the lower bound becomes unity instead of zero. This is necessary since CO = 0. If CO + 0, then the sum would become infinite at n = 0. A necessary requirement for all the coefficients is that they be finite. Therefore, the recursion formula is: for n > 0 5C C 2 (0) 6-ukCn (2)' C =(I - 5k - ak2) C (0)+an_ C + 2ak n (40) ~n n na n,O na where 0o Ck C (2) ~~ Ck k k-0 a k(n - k) k= 0 For n = 0, we determine k as follows ck2 + 6k- =O (41) or k=- (42) For n = 1, 5C C 1 (2) C ac C + 2 (ak 1 a 1,0 a or C (-+. =0 (43) The second factor is not zero; hence C1 =0 (44) For n = 2, it is possible to determine the characteristic exponent a: 2 6C2 aC1 2 akC 2 1 2 (45) C 1 - (45) 2 2a 2 2a a - =0 (46) \a a 2a/ 13

Institute of Science and Technology The University of Michigan Either or both factors will be zero. We shall assume that the second factor is zero to determine a. By taking C1 = 0 from Equation 46 and substituting the value of k given by Equation 41, we get the result: a = ~ 1/2 V56 + 4ao (47) This result is known to be correct since it can be compared with the solution obtained by separating the variables (see Appendix A). The characteristic exponent here has two values. Both must be used in Equation 36. This gives rise to two Dirichlet series. If the characteristic exponent had k values, then k series would be needed in addition to the, constant of integration which may arise. This procedure is similar to the one used in linear differential equations. The remainder of the solution follows a procedure similar to that described elsewhere [13]. If the second partition scheme is used, a different form of x(t) is obtained; this is easily checked. Also, the equation of the auxiliary forcing term for A(t) in the forms A(t) =' Cn tneat (48) n=O A(t) = Cn tnenat(49) n=O =0 A(t)= C tn (50) /A! n n=O eventually follows the same procedure and yields different and interesting cases for investigations. These cases enable the engineer to gain further insight into the nonlinear system under consideration. Suppose, e.g., the form of Equation 50 is chosen and the following results are obtained for the same system described by Equation 32. The moments are tn+1 Qn(t)= n~1 (51) Thus G C x(t) = ntn+ (52) n=O and it is possible to obtain the recurrence relations C C (0) n C (2) (53) n n n n-2 14

Institute of Science and Technology The University of Michigan Hence the first few coefficients are CO = 3 C1 = -6j (54) C2 2 3 46~2 C A9 + 4ba3 3 1 3 Professor S. S. L. Chang [42] (New York University) suggested that the partition should be made so that N, the nonlinear term, is the smallest possible. Professor Chang offered the following illustration: 0.01 xR + sinh x + x = 0 (55) For this equation, it is better to use Z(D) = 0.01 D 2+D+1 (56) N = (sinh x - k) than Z(D) = 0.01 D + 1 (57) N = sinh k It is important to mention here that the types of nonlinearities which can be solved by use of the partition theory can be divided into two classes: (1) Algebraic Nonlinearities: These are exemplified by a polynomial in x and its admissible derivatives, and they are described by: N N(x) = akx (58) k=0 M N m N(x)= L Z amk{X( (59) m=O k=0 M N X= aakx(m)x(k) (60) m=0 k=0 where x = x(t) and the superscript (k) indicates the k-th derivative of x with respect to time. 15

Institute of Science and Technology The University of Michigan (2) Transcendental Nonlinearities: These are restricted to functions of x(t) which have an inverse whose derivative is a rational function, so that the system equation can be easily transformed to equations with algebraic nonlinearities. Typical examples are Z(D)x(t) + bex(t) = g(t) (61) Z(D)x(t) + C tan h x(t) = g(t) (62) J 0'+ r 0+ a( sin 0= T (63) where b, c, J, r, (a, and T are constants. Equation 61 occurs in problems dealing with physical electronics, Equation 62 in physical systems containing saturating types of nonlinearities, and Equation 63 in synchronous machinery. To show how these types of equations can be transformed into equations containing algebraic nonlinearities, let us take Equation 61 as an example. If we assume W = ex (64) then x = log W (65) Substituting Equations 65 and 64 into Equation 61 gives Z(D) log (W) + bW = g(t) (66) Differentiating Equation 66 with respect to (t) gives Z(D) W + bW = (t) (67) if we assume that Z(D) = Y(D) (68) Substituting Equation 68 into Equation 67 gives w - W[Y(D)g(t)] + bW [ Y(d)x*J] = 0 (69) Certainly, Equation 69 cannot be solved for all forms of Z(D) by the method presented, because the term Y(D)WV which appears in the equation is difficult to handle. Moreover, Equation 69, as given, is an integral differential equation. For the special case in which Z(D) is the firstorder differential operator, Equation 69 reduces to the following form: W - 1[*(t)] w + bW D() = O (70) 16

Institute of Science and Technology The University of Michigan When the operator 1/D is replaced by an integration operation, Equation 70 yields W - g(t)W + bW2 = 0 (71) Equation 71 can be solved by partitioning techniques. When Z(D) is a linear combination of x and one or more of its derivatives, the problem is more complicated but can be explained by the case in which dx Z(D)x +x (72) Substituting for Z(D) from Equation 72 into Equation 67 gives (D + 1) + bW = g(t) (73) or d-t (-+ W + bW = g(t) (74) +W + bW = g(t) (75) W2 Equation 75 can be put in the form: W'+ - (t)W + bWW = 0 (76) Equation 76 is of the form: Z(D)W + N(W, W) = 0 (77) where Z(D)W= [D2 +D-(t)]W and N(W, ) =bW2W -W2/W (78) Further investigations of N(W, WIV) show that it satisfies the Lipschitz condition. A review of A. A. Wolf's published work, by Professor L. F. Kazda and the author [43], revealed that certain questions have been left unanswered, which the novice, when trying to apply the method presented, would naturally like to have summarized. These questions are: (1) What insight can be gained about the stability of a system for a general class of forcing functions by using the partitioning method of analysis? 17

Institute of Science and Technology The University of Michigan (2) Under what conditions would partitioning at points other than the highest derivative be desirable? (3) How does the partitioning method help one to understand the behavior of the amplitude and frequency of oscillations with time? Consider, e.g., the following nonlinear differential equation: x + wx + x = 0 (79) Reference 43 gives answers to the above questions. In the author's opinion the partitioning technique has not been used in the field of nonlinear synthesis as much as in the field of nonlinear analysis. However, Professor Y. H. Ku (Moore School of Electrical Engineering, University of Pennsylvania).has pointed out that it is precisely in the field of synthesis rather than the field of analysis that the new method is most useful. Dr. A. A. Wolf [43] has suggested a synthesis procedure as follows. Let Y(s) be the transfer function of the linear part and let {sn be a set of singularities of the system response in the complex plane generally supposed to be finite in number. Then, given a specified forcing function g(t), what is the nature of the nonlinear part N that enables it to satisfy these conditions? The synthesis procedure would be as follows: (1) From the singularities and initial conditions, determine C in closed form as shown in Reference 2. (2) Use the C coefficients to determine the response of the system by using the moment theorem. (3) Then determine the nonlinear part of the system by using the canonical Equation 1. (4) Finally, by considering Figure 6, the system is obtained easily. As an example [43], let the linear part of a system be given as: Y(s) = s (80) If the system is forced by g(t) = e-2t (81) find the coefficients of the nonlinear terms so that the response has a pole at s = -1 and the nonlinear part has no derivatives and does not exceed the second degree. The solution can be 18

Institute of Science and Technology The University of Michigan effected by using the recurrence equations of References 45, 46, and 47, or it may be solved directly. In the latter case, the solution is almost trivial. From Z(s) Y(s) + 2 (82) the required system equation form is dx/dt = 2x + N(x) = g(t) (83) N(x) = g(t) - i - 2x (84) Note that N(x)= akx (85) k=O so that E akx = g(t) - x - 2x (86) k=O From the problem, it is given that X(s)= 1 (87) from which x(t) = e-t (88) so that -t + a 2t -2t -t (89) a +ad +ae- (89) From the identity of Equation 89 it is clear that a0 = 0 (90) a= -1 (91) a2=1 (92) Thus the nonlinear component satisfying the problem is N(x) = -x + x2 (93) 19

Institute of Science and Technology The University of Michigan 1.2. SUMMARY AND RESEARCH OBJECTIVES A short review of a generalized mathematical theory for analyzing a certain class of nonlinear systems has been presented. It was stated that under a broad set of conditions imposed on the linear, nonlinear, and forcing functions, a class of systems exists for which it is possible to develop exact and unique mathematical solutions. Although linear differential equations generally have unique solutions, examples show that nonlinear differential equations generally do not have such properties. The techniques of partitioning for systematic solution of actual physical problems in connection with feedback systems also were reviewed. Examples were shown to illustrate that a power series solution of differential equations is just a special case of the partitioning technique. It is interesting to note that a complete synthesis theory can be developed by using the system shown in Figure 6 as a basis. In addition, it is possible to show that by appropriate partitioning almost any desired configuration can be obtained in both single and multiloop feedback systems, a fact of paramount importance in engineering today. For these reasons, the author felt that it would be of interest to control engineers to extend WolfIs work to include this new class of nonlinear systems. The objectives of this extension are (1) To develop, by using the partitioning technique, a theory for the analysis of a certain class of systems, the behavior of which can be described by a single-degree-of-freedom nonlinear differential equation, with time-varying parameters. (2) To apply the partitioning technique and the state-variable approach in further investigation of the system's behavior and stability. (3) To extend items (1) and (2) above to include systems with more than one degree of freedom. (4) To analyze the system numerically by using a digital computer for systematic approximations of previous theoretical results. The Runge-Kutta methods of analysis are to be adopted for the purpose of comparison. 2 A THEORY FOR THE ANALYSIS OF A CERTAIN CLASS OF NONLINEAR TIME-VARYING PHYSICAL SYSTEMS The partitioning technique introduced in Section 1 is used here to analyze a certain class of nonlinear physical systems, whose dynamic behavior can be represented by a nonlinear differential equation having some time-varying terms and a forcing function. Most of the material in this section describes work done under the guidance of Professor L. F. Kazda while the 20

Institute of Science and Technology The University of Michigan author was at The University of Michigan [2]. It will be shown that by placing certain restrictions on the system equations a unique solution which belongs to L2 space can be obtained.5 This result is useful because the analysis presented gives a means of controlling the amount of energy associated with the system response; but, it is also limiting because the shape of an L2 function may not be desirable in all applications. (This disadvantage will be discussed in Section 3.) The exact analytic solution is given as the limit of a sequence of Picard iterates {Xn} each of which is also an L2 function. Uniform convergence of these iterates is shown to be guaranteed by the above restrictions. The well-known mean square error criterion is then applied to find the number of iterates necessary to approximate the solution of this class of nonlinear differential equations. An upper-bound function P(t), which also belongs to L2 space, is calculated. Variations of the solution with respect to initial conditions and system parameters are investigated. Finally, some applications to specific nonlinear problems will be displayed to illustrate the method presented. 2.1, SYSTEM DESCRIPTION AND DEFINITIONS Consider a certain class of physical systems that can be described by the general, nonlinear, time-varying differential equation of order n, namely: G(x,,,..., x(n)t)= 0 (94) where x(t) and its first (n - 1) derivatives are given att = 0. Furthermore, assume that Equation 94 can be put in the form [am(t)Dm] x(t) + N(x,k...., t) = g(t) (95) m=O where [D] = d/dt = the differential operator am(t), s = time-varying parameters g(t) = the forcing function6 N(x, xi,..., t) = the nonlinear term t = the independent variable Equation 95 can be written in the following form: L(D, t)x(t) + N(x) = g(t) (96) 6For an explanation of L2 spaces, see References 4, 30, 38, and 40. Also, Appendix B includes information on some of the properties of L2 spaces which are considered in this report. 6 See Appendix C for situations in which there are forcing functions of the form K(D, t)g(t). 21

Institute of Science and Technology The University of Michigan where L(D, t) is the linear operator operating on the system output x(t) and is given by L(D, t)x(t) = Z [am(t)Dm]x(t) (97) m=0 The following definitions [17, 48, 49] can now be introduced: (1) Let W(t, u) represent the impulsive response of the linear part of Equation 96 which is defined by Equation 98:7 L(D, t) W(t, u) = 6(t - u) (98) where 5(t - u) is the unit impulse applied at t = u. (2) Let x0(t) be the response of the linear part with g(t) as the input forcing function and with the same initial conditions as in Equation 94. If zero initial conditions are assumed, x0(t) is defined as follows: t x0(t) = { W(t, u) g(u)du (99) 0 This means that Equation 99 represents only the particular solution of the linear equation L(D, t)x(t) = g(t). In general, the analysis is true when x0(t) is the complete solution which contains the effect of initial conditions of the original nonlinear system and when the system restrictions given in Section 2.2 are satisfied. (3) Let F[t, u, x(t)] = W(t, u) N(x,., x(-), t) (100) be a nonlinear function obtained by combining the actual system nonlinearity N and the weighting function W(t, u); hence, it can be called a weighted nonlinear function. In general, F depends on x, x,..., and x(n-l) 2.2. MATHEMATICAL RESTRICTIONS AND THE SYSTEM ANALYSIS In order to present a theory that applies to a general class of problems, a number of mathematical restrictions are now placed on the system parameters and on x(t) and F(t, u, x), which were defined in Section 2.1. For simplicity of analysis, the weighted nonlinear function F is assumed to be a function of t, u, and x and is not.dependent on c, x,.... For example, F can assume the following forms:8 7See also Appendix C, 8The extension to more complex nonlinearities is explained in Appendix F. 22

Institute of Science and Technology The University of Michigan F(t, u, x) W(t, u) w(t,k k=l (101) F(t, u, x) - W(t, u) [exp (x)] The restrictions required are: (1) xO(t) is an L2 function9 in the interval of interest [0, T.] 0 < T < + coo (102) (2) For a given region D, [Ix I < d(t), 0 < u < t < T], the function F(t, u, x) satisfies the following conditions: (a) F(t, u, x) satisfies the following conditions; namely, if the two triples (t, u, Zl) and (t, u, z2) are in the domain D, then IF(t, u, z1) - F(t, u, z2)I < K(t, u) z1 - z21 (103) where K(t, u) is an L2 function in the given domain, that is, 0 K (t, u) du < a2(t) 0 < t < T 0 and (104) T 2 a 2(t) dt < A = constant 0 (b) j IF[t, u xO(u)] du < n(t) (105) where n(t) belongs to L2[0, T], that is, n2 (t) dt< C = con stant (106) 0 and A2 and C2 are two finite positive numbers. (c) It is required that: x0O(t) < d(t) (107) q(t) < d(t) 9The definition of an L2 function on the real line and on the plane is given in Appendix Bo Also, restriction (1) can be written as: x0(t) belongs to L210, T]. 23

Institute of Science and Technology The Ciniversity of Michigan where t x1(t) - x 0(t) F(t, u, x0(u)) du and (108) cc k oo Ak q(t) - xl1(t)l + C c(t)E k=0 where A, C, and ac(t) are as defined before. (3) The L2 functions x0(t), n(t), and oa(t); the time-varying parameters ai(t); and the enforced bound d(t) are assumed continuous and bounded over [0, T]. (109) It should be noted that the function d(t), appearing in the definition of the domain D, may or may not be an L2 function. This depends on the particular problem. Its use as a function of t gives flexibility to an analysis of problems in which the interval of interest is infinite. Unfortunately, it is not helpful for all problems, as illustrated by the following special examples: Let F(t, u, x) = W(t, u) N(x) = eX(t)x (110) Then IF(t, u, z1) - F(t, u, z2)1 = e-(tu)z1- z2< K(t, u)z1 - z 2 (111) Therefore, K(t, u) > eX(tu) (112) Thus, in this particular case, K(t, u) is independent of the choice of d(t) and does not belong to L2 space for infinite'intervals, but only for finite intervals. Consequently, for problems involving infinite intervals, the linear terms of x should be included in the linear operator L(D, t)x, and not in N(x). Again, assuming that cc N(x) = anx (113) n=2 then F(t, u, z1) - F(t, u, z2) = W(t, u)(z1 z2) an[zln + 2 2 n (114) n=2 24

Institute of Science and Technology The University of Michigan Letting I z < d(t) > 0 gives IF(t, u, z1) - F(t, u, z2) = IW(t, u) | an n Z 1 + + z2 z1 + 2 +z Z1 z2 n=2 (115) Defining: Q(t) =max a n-1 + n-2 + + z n-2 + n-1 (116) d t) ~+ a 1 z2z1 +. 2 1 + z2 z1 5. d(t) n=2 z2 < d(t) gives IF(t, u, z1) - F(t, u, z2)I < W(t, u)Q(u) z1 - z2 = K(t, u)z1 z2 (117) Hence, K(t, u) = W(t, u)Q(u) (118) Therefore, for problems with infinite intervals, d(t) and hence Q(t) should be chosen so that LJ K (t, u)du] = w (t, u) Q (u)d1 = c(t) e L2[0, L 0] (119) Now, applying the partitioning technique to the system described in Equation 96 yields L(D, t) x(t) = g(t) - N[x(t)] = f(t) (120) Equation 120 can be considered as a linear differential equation with f(t) as an auxiliary forcing function. This function is not unique since it depends on how Equation 96 is partitioned. For a particular choice of partition, N(x) may include some of the linear terms of Equation 96.O Equation 120 is known as the linear partitioned system. If zero initial conditions are assumed, x(t) can be obtained from t x(t) = J W(t, u) f(u) du (121)' See the rules of partitioning in Section 1, or see Reference 47. 25

Institute of Science and Technology The University of Michigan or, equivalently, from t x(t): = W(t, u) {g(u) - N[x(u)]}du 0 t x(t) W(t, u) g(u)du - W(t, u) N[x(u)]du (122) 0 x(t) = x0(t) - $ F[t, u, x(u)]du Equations 96 and 122 are equivalent, and x0(t) and F are as previously defined. Also, Equation 122 still holds for problems with nonzero initial conditions when x0(t) is the complete solution of the linear equation L(D, t)x0(t) = g(t), which contains the effect of these initial conditions. Equation 122 is a nonlinear integral equation of the Volterra type [26, 38, 39], and it can be shown that for the class of systems represented by Equation 96 and the given system restrictions, theorems A, B, C, and D (which follow) are true. In general, the proof of theorems A and C follows the method of Vito Volterra for dealing with nonlinear integral equations [38]. Since the proofs are important for an understanding of the next sections, they are given in detail. 2.3. THEOREM A Let a set of functions x0, x1,..., xn, each of which is a functionof time, be defined according to the following recurrence relation: xn+ (t) = x0(t) - F[t, u, xn(u)]du (123) where n = 0, 1, 2,..., and x(t) is defined as in Equation 99. Then, under the restrictions placed on the system (Equations 102 through 109), the sequence {Xn} belongs to an L2 space and is in the domain D -[Ixl < d(t), 0 < u < t < T]. Proof It can be shown that when t lies in the interval [0, T], the sequence {Xn} stays inside D. Assume (as an induction hypothesis) that xi(t)(t), t),.., Xn(t) are inside D; then Equation 123 shows that x n+(t) is defined. Replacing n + 1 by n in Equation 123 gives Xn(t) = x0o(t) - F[t, u, xnl(u)]du (124) 26

Institute of Science and Technology The University of Michigan Subtracting corresponding sides of Equations 123 and 124 and squaring gives [Xn+(t) - x n(t)] = F[t, u, xn(u)] du -j F[t, u, x 1(u)]du = F[t, u, x (u)] - F[t, u, x 1(u)] I du] (125) Making use of the restrictions gives t 2 [Xn+(t) - x (t)]2 < i K(t, u) Ixn(u) - x 1(u) Idu (126) Applying Schwarz's inequality to Equation 126 gives [Xn+l2(t) Xn(t) < K(t, u)du X n(u) - x (u) du < o(t) t[xn(u) - xn l(u)] du (127) J0 for n = 1, 2,.... For n =, Equation 127 gives [x2(t) - x1(t)]2 < a2(t) (x1 - x0o)2 du (128) The integral on the right-hand side of Equation 128 can be obtained from Equation 123 by putting n = 0, as follows: x1(t) = x0(t) - F[t, u, x0(u)]du (129) 0 Equation 129 gives T - 2 fT n2(t)dt 2(130) 0 [xl(t) - xO(t)dt< (t) dt < c (130) Thus, x2(t) x1(t)]2 C2 2 [x2(t) - x(t)] < C (t) (131) 27

Institute of Science and Technology The University of Michigan and (x -x) <C C2A2 (132) 0 (x2 -x1< Appendix D shows that 2 2n [xn (t ) x(t)]2 dt (133) Equation 133 gives Ix (t) - x (t)]2 dt < C2 A(n ) (134) 0 n n-i dt - (n- 1)! Substituting Equation 134 into Equation 127 gives 2 2 2 A( [x n+1 l(t)- Xn(t) < C (t) (n -1)' (135) for n = 1, 2,... Now, ixn+l(t) = lxl(t) + [x2(t) - x1(t)] + [x3(t) - x2(t)] +... +.. + [xn+1(t) - Xn(t)] < Ixl(t) + Ix2(t) - x1(t) l + * *. + I t)n+ (t) (136) Substituting from Equation 135 into Equation 136 gives [x + (t) I < I x(t) I + C A(t) + C ( IXn+l(t)1 < LXl(t + c c27(t) c cn(t) + c c~(t)-~,. +... +... + c c(t) ~n _ ii) 00 Ak (137) = 1x1(t) I + C a(t) = q(t) (137) k=O But, condition 108 shows that q(t) < d(t); thus, Ix (t) < d(t)0 < t < T (138) Since Ix (t) and x1(t) are in D, the sequence {xn(tj} given by Equation 123 is meaningful and is also in D. Note that x0(t) is assumed to belong to L2[O, T] and Equation 130 shows that x1(t) also belongs to L2[0, T]; hence Equation 135 shows that the sequence {Xn(t)} is in L2[0, T] for n = 1, 2, 3,... 28

Institute of Science and Technology The University of Michigan 2.4. THEOREM B If conditions 102 through 109 are satisfied, then the sequence of functions x0, x, x2,... x of Theorem A converges in the quadratic mean to a function qp(t).It also converges uniformly in the ordinary sense to a function t4(t). The two functions Q((t) and g/(t) are identically equal to a function x(t) for all t in the region 0 < t < T, where T is the period of interest. Proof For convergence in the mean, we have to prove first that, given two numbers m and n so that m > n, oJLim [Xm(u) - xn(u)] du = 0 (139) m-oon Now, since (xm n ( x - x + ( + + (x + + -) (140) -x)=(x Xnmm-1 m-1 -2 n+1 n therefore k=n Ix - x I < IXk kl (141) k=m-1 Substituting Equation 134 into Equation 141 gives k=n (k-l) Ix - x I < C a(t) (142) m n(k - 1 k=m -1 Substituting Equation 142 into the left-hand side of Equation 139 gives T T k=n K-12d Lim X [xm(u) -X(u)] du = Lim CJ k du n-co m n- n-o 0.- m —o 1)2 im C22 (k-1) Lim C 2A =m0 A (143) n- -1) -k=m-1 Therefore, Equation 139 is satisfied. Since xm is an L2 function, x is also an L2 function; hence, by using the first theorem given in Appendix B one can conclude that there exists a function p(t) of L2 space to which the sequence {Xn } converges in the mean, and this function op(t) satisfies the following equation: Lim i [xn(U) - (u)] du = 0 (144) 29

Institute of Science and Technology The University of Michigan or, in shortened notation, Lim x (u) = sp(u) (145) To prove convergence in the ordinary- sense, one can form the series: x1 + (x2 -x1) + (x3 -2) + (146) In this series, n-2 ( Xn x 1) < C a(t) (n A 2) (147) Since a(t) and x1(t) are assumed finite, one can easily conclude the required uniform convergence by the M-test. But the n-th partial sum of series 146 is given by =x + (x2 -x1)+ + (x-x ) (148) Therefore as n - co, x has a limit. n Let Lim x (t) = (t) (149) n-Hc n Note that uniform convergence of the sequence xn can be easily concluded from Equation 142, since this equation shows that the sequence x is a Cauchy sequence in the ordinary sense. n Now, by applying H. Weyl's lemma (stated in Appendix B), we can conclude that 4(t) = y(t) (150) almost everywhere.'Let P(t) = co(t) = x(t), 0 < t < T (151) Equation 150 holds everywhere because the functions are continuous. The existence, uniform continuity, and boundedness of the first n derivatives of x(t) are discussed in theorem E of Section 3. 2.5. THEOREM C Under the conditions placed on the linear, nonlinear, and forcing functions of the system given by Equation 96, the function x(t) as defined in Equation 151 is a solution of the differential Equation 96, and is the unique solution for a given initial state. 30

Institute of Science and Technology The University of Michigan Proof To prove that x(t) is a solution of the differential Equation 96, one can easily form the following equation: [Xn - (t)] (t) ) x0(t) + F[t, u, x l(U)]du (152) Integrating both sides of Equation 152 in the interval (0, T) and using Equation 144 gives Lim T [X- (t)]2 dt = Lim fT [(t) - x0(t)] + F[t, u, x (U)]d = (153) Since L.i.m. x (t) =L.Ai.m. x (t) = (t) (145) n —:c& n-i n-cc n therefore, O(t) - x0(t) + F[t, u, q'(u)]du = 0 (154) hence, p(t) satisfies the integral Equation 125. Consequently, x(t) is a solution of the differential Equation 96. To prove uniqueness, assume (o(t)* to be another solution for Equation 96. Therefore, (t) (t) xt + F[t, u, sp(u)]du (155) o* (t) = x0(t) + F[t, u, o (u)]du (156) Subtracting Equations 156 and 155 and squaring gives [o* (t)- (p(t)]2 = ( {F[t, u, (t)] - F[t, u, 9 (u)}du (157) Applying the restrictions on the system gives K It T * (t) - o(t)] < K (t, u) du [*(u)] du < a (t) [*(u) - ~(u)] du (158) 0 0uu Let [s*(t) - o(t)]2dt = S (159) 31

Institute of Science and Technology The University of Michigan Therefore, [(*(t) - c(t)] <S a (t) (160) Substituting Equation 160 into Equation 158 and repeating n timesll gives T 2 2 A2n(*) [p*(t)- (t) dt < S n' (161) Letting n - oo gives T] {T [*(t) 0(t2dt = 0 (162) 0 Hence eo = CD, and the solution x(t) of the system equation is unique. 2.6. THE APPROXIMATION OF THE SOLUTION IN THE MEAN SQUARE SENSE According to the preceding theorems, the solution x(t) of the system Equation 96 was found to be the limit of a sequence of iterations {Xn} in the L2 space as the value of n tends to infinity. It may be reasonable to ask whether there exists a number H such that for all n > H, the error Er(t) remains within a specified value. It is apparent from the mathematical definition that the number H exists. A reasonable error criterion for the calculation of H is the mean square error criterion. Thus, one can form the following equation: n2 (u) d(u) < q, for n > H (163) where 2 = (x - x )2 (164) 2 xn is an approximation for the actual solution x(t) of the system under consideration, and q is a prescribed positive error. It should be noted that the analysis also leads to an estimate of En itself. The author uses the mean square error in the analysis, because this leads to a result which is simple to use in the applications. "After the first substitution, [o* - ]2 S2a(t) 2(t) 2(z)dz. By repeated substitution, [v*(t) - sp(t)]2 < S2a2(t) | a2(y)dY a2(z) dz.... Integrating both sides of this equation gives Equation 161. 32

Institute of Science and Technology The University of Michigan Now, by making use of the following equation: k=n (k 1) lx -x < Ca(t)A (142) k=m-1 holding n fixed, and letting m - co, one gets oo k-i n= ix - xnI < E (tA- (165) k=n Therefore, _ < Ca(t) ( -A) (166) n - k=n-1 2n EC ( (167) =n-1 Substituting Equation 167 into Equation 163 and making use of Equation 104 gives: r n (u)=du < C A q (168) J )0 CA (k-d 2 However, Equation 168 is very difficult to solve for n if constants C, A2, and q are specified. But, a simplification for this situation can be obtained by making use of the following inequality [8, p. 13]: 0o Ak An- 1 oA k,A.< A E. (169) k=n-1 -k=0 Substituting inequality 169 into Equation 168 gives 2 2 2n 2 C ( O A (170) k=0 Equation 170 gives an upper bound for the error when the actual solution x(t) is approximated by the function Xn, in the sense of the mean square. The infinite series of Equation 170 can be shown to converge for all finite values of A. In addition, as the number of iterates n tends to infinity, the square of the error q tends to zero. This condition is explained by theorem C since 33

Institute of Science and Technology The Oniversity of Michigan the actual accurate solution can then be obtained. In general, the resulting square error q2 depends on the number of iterates n and on the constants A2 and C2, already defined by Equations 2 104 and 106 respectively. Therefore, for a given mean square error q, the number of iterates 2 2 necessary to approximate the solution is dependent upon A and C and, hence, upon how the system Equation 96 is partitioned. When Picard iterates are used to yield the approximate solution of the system, it is desirable to partition the equation so that fewer iterates are needed for a given mean square error q. This reduces the labor required for the solution of the system equation. Unfortunately, it is very difficult to formulate a rule for choosing the partitioning so as to obtain the minimum number of iterates. However, for a given n-th order equation, there are an infinite number of choices for partitioning the equation.12 At each choice, the number of iterates (n) necessary to yield a mean square error q is obtained by using Equation 170 after calculating A and C. In order to compare two partitionings, let us assume that for the same q1 we have C1, A, nI and C2, A2, n2 as the C, A, n numbers for the first and second partitionings, respectively. Applying -Equation 170 yields 2n 2 2n2 2 2(n2 - 1 A 1A27 A C ( 1 ( 1 -1) A22n2- 1)'. 2 kk= =O k or 2n1 (%n- 1)' A1 P - 1 2 (172) (nI P. ~2n2 Pi where k 2 P k=C and (173) P2 = C22. 12See the rules of partitioning, Section 1. 34

Institute of Science and Technology The University of Michigan 2 An examination of Equation 172 shows that in order to have n2 > nl, A2 should be greater than A and C2 should be greater than C1. 1 2 1 When A is less then unity, the infinite series of Equation 170 will be quickly convergent. 2 2 2 For a given A2, C2, and q, the value of n that satisfies Equation 170 is best obtained by the following graphical procedure. Consideration of Equation 170 gives 2n 2 A q (174) (n - ) 2 k=0 or 2n A2n - Q 2 (175) (n- 1)I where Q2= q 2 (176) c2 Ak) k=0 ~,.L-11X-IU A V1 L-LYLLLJIL 174 can be drawn for different values of n by using a given value of A. A straight line drawn parallel to the n axis from a distance Q intersects the last drawn curve at the required point. As an illustration, the left-hand side of Equation 174 is drawn for the cases in which A = 1 and A = 0.364 < 1 in Figure 7. The case in which A = 2 > 1 is also plotted in Figure 7 but, for convenience, not to the same scale. The Figure shows that if the values of Q 2 2 A2n as defined by Equation 176, are chosen so that Q < A = (n- 1) n=l' then there will be only one value of n corresponding to a given mean square value of error. This follows because a straight l r i 2 A2 line drawn parallel to the n axis at a distance Qfrom it will intersect the curve (n - 1) at one 2 2 2 A2n and only one point. If Q is chosen so that Q > A = (n - 1) then the magnitude of A2 will determine which of the following situations will occur: 2n (1) If A2 < 1, then no solution will give a positive value of n because the curve n is ~2 ~~(n-i monotonic decreasing, and has a maximum value which is equal to A and occurs at n = 1. Thus a straight line drawn parallel to the n axis with Q > A will not intersect the curve (n for any positive value of n. 35

titute of Science and Technology The University of Michigan 2 A2n 1.6- / —— A = 2.4 - -A=.364.2 - 0 2 3 4 5 6 7 8 9 10 11 n 2n FIGURE 7. A PLOT OF VERSUS n (n - 1); 2n (2) If A2> 1, then the corresponding curve (n - 1) increases as n increases, until it attains its maximum value at a particular value of n. With this behavior of the function 2n A, one expects the following two situations to happen: (n - 1)' 2 A2n (a) If Q > max (n - 1) then no solution will give a particular value of n, because A2n 2 there will be no intersection of the curve with the straight line Q drawn parallel (n - 1)! to the n axis. ~~~~2,~A2n (b) If Q2 lies between the value of the function An at n = 1 and its maximum value, 2n 2 2 A 2 that is A< Q < max then there will be two values of n for a given value of Q n (n-) that is, for a given mean square error. This is not surprising because it indicates that among the sequence (Xn}, there are two elements xnl and xn2 (n1 < n2) which will approximate the solution x(t) for the same mean square error. For the purpose of simplicity, one should select xn1 which has the smaller index n1. 2.7. CALCULATION OF AN UPPER-BOUND FUNCTION:(t) FOR THE SYSTEM RESPONSE x(t) To calculate an upper bound:(t) for the system response x(t), one should note that the sequence {xn(t)} is uniformly bounded by q(t). This is given by the following equation: o k Ix l(t)i (t)l < q x(t) = ix (t) + C (t) (137) 36

Institute of Science and Technology The University of Michigan for n = 1, 2,.., co. But theorem B shows that Lim x (t) = x(t) (177) n-co n uniformly. Therefore, lx(t) i < q(t) (178) Also, since xl(t) = x0(t) - { F[t, u, x0(u)]du (179) it follows that xl(t)l < 1x9(t) l + In(t)I (180) Now, define a function p(t) so that p(t) = [x0(t)I + In(t) I + mac(t) (181) where co k m = c (182) 0 and compare Equations 137, 180, and 181 which give q(t) <:(t) (183) The upper-bound function:(t) should obviously stay inside the domain bounded by the enforced bound d(t).'3 0 < lx(t)I < q(t) < P(t) < d(t) (184) From Equation 184 it is obvious that both q(t) and:(t) work as upper bounds for the system response x(t) inside the enforced bound d(t). Although q(t) is a lower bound than P(t), and is a better bound from this point of view, the latter has the advantage that is may easily be obtained in calculations, as shown in Equation 181. 2.8. VARIATION OF THE SOLUTION WITH RESPECT TO THE INITIAL CONDITIONS AND THE SYSTEM PARAMETERS (THEOREM D) Let the initial conditions x(i)(o0) change to x(i)(o) + ( tx(i(), ( i = 0, 1,..., n - 1), then x0(t), which was defined by Equation 99 and which in addition contains the effect of the initial 13See Figure 12. 37

Institute of Science and Technology The University of Michigan conditions, will change to [xO(t) + w(t)]. Instead of the solution x(t), there will arise a new solution X(t), as given by X(t) = [x0(t) + w(t)] - F[t, u, X(u)]du (185) If we require that Ax(i)(0) = t so that 11 < 6, and that this change in the initial conditions will keep [w(t)] e L2 in the interval of interest,14 then the existence.and uniqueness of the system Equation 185 follows exactly as before. This particular change of initial conditions of the original nonlinear system requires an investigation of the behavior of the norm of the difference of the two solutions X(t) and x(t), symbolized as II x(t) - x(t) 11.' This investigation is carried on in the following theorem. Theorem D Under the conditions (given by Equations 102 through 109 placed on the linear, nonlinear, and forcing function terms of the system described in Equation 96, the system solution is uniformly continuous with respect to the initial conditions x(i)(o). If the system contains a parameter h, then the solution also is uniformly continuous with respect to X. Therefore, x(t, x(i)(0), X) is uniformly continuous with respect to x(i)(0) and X. Proof Subtracting Equations 185 and 124 gives S {*t [X(t) - x(t)] = w(t) -IJ (F[t, u, X(u)] 1 F[t, u, x(u)] du1 (186) Applying the Lipschitz condition to Equation 186 gives [X(t) - x(t)] = [w(t)] - {Kt(u) [X(u) - x(u)]}du (187) 0 Squaring both sides and using the Schwarz inequality gives 2 [X(t) - x(t)] < 2{w(t) + a (t) j [X(u) - x(u)]du (188) 4Since x0(t) e L2[O, T], and w(t) e L2[0, T]; this implies that [x0(t) + w(t)] e L2[0, T]. =X(t) x(t)II2 J(X - x)2 dt 0 38

Institute of Science and Technology The University of Michigan Letting S2 t 2 Letting S [X(u) - x(u)] du, and substituting into Equation 188 gives [X(t) - x(t)]2 < 2[w(t)2 + S2a 2(t)] (189) Substituting Equation 189 into Equation 188 and repeating the first substitution gives [X(t)- x(t)] 2 2 w(t) + 2 2 (1 + 2 + 22 (2)2 23 2(3A 2n322(tn +2 222n (A2n + + 2. )j + 2a(t)S n (190) where 2 St w2(t) dt < oo 0 The integral appearing in Equation 191 is finite since w(t) E L2 by definition, 2 can ar be made small by a proper choice of Ax(i)(0). As n tends to infinity the last term in El 191 drops out and we are left with [X(t) - x(t)] < 2 L (t) + 2 2 2 I (2A) n=0 or [X(t) - x(t)]2 < 2 [w2(t) + 22a2e2A e (192) Integrating both sides of Equation 192 in the interval 0 < t < T gives IIX(t) - x(t) 12 < 22 [1 + (2A2) e(2A2 (193) or I1X(t) - x(t) 1]2 < 2 (194) where 2 = 22 [1 + (2A2) e(2A2] (195) Equation 194 shows that the two solutions x(t) and X(t) may be made to differ by a small quantity, as we please. A proof proceeding on similar lines as the above shows that if the system equation involves a parameter X, then the solution x(t) is continuous with respect to X. Hence the theorem is true. 39

Institute of Science and Technology The University of Michigan 2.9. THREE EXAMPLE CALCULATIONS Example One Let us consider a system whose behavior is governed by the following nonlinear differential equation: dx 2 -2t dx + x + x = e (t = 0, x = -1) (196) This can represent either the electric circuit of Figure 8 or the nonlinear feedback systems of Figure 9. Now, partitioning Equation 196 at the second term gives dx -2t 2 dt + x = e2t x2 = f(t) (197) The impulse response of the left-hand side of Equation 197 is obtained by solving the following differential equation: dx dt + x = 6(t - u) (198) This solution is W(t, u) = W(t, u) = +e(tu) (199) The term x0, as defined before, is the solution of the linear equation: dx -2t + x = e (t = 0, x=-1) (200) Therefore, x0(t) = -e2t Now, since F[t, u, x(t)] = e (t-u)x2 (201) R=1 L=1 e-2t R(i) = i FIGURE 8. AN RL CIRCUIT WITH A NONLINEAR RESISTOR 40

Institute of Science and Technology The University of Michigan ~~~~~~e + 1~~ x(t) s+l x 2 x(t) N(x) (a) e +,y x(t) (b) e-2t +~ "-{ x(t) =e ~~~~1/S 2 l (c) FIGURE 9. EQUIVALENT REPRESENTATIONS OF THE SYSTEM GOVERNEP BY EQUATION 185 41

Institute of Science and Technology The University of Michigan then F[t, u, x0(t)] =e (tU)x 2 e-(t-u)-4t (202) Applying condition 105 of the restrictions to Equation 202 gives n(t) n(t) = F[t, u, x0(u)]du = j e(U)e -4u 1/3 (e -e 4t (203) Also F(t, u, x) satisfies the Lipschitz condition in the domain D = [x < 1, (0 < u < t < T], since F(t, u, x1) - F(t, u, x2) = e(tu)x1 - x2] = e (t) x1 - x2 I1 + x21 < 2x e (txu)x - x (204) Hence IF(t, u, x1) - F(t, u, x2)1 < 2e(tU)x - x2 1 (205) Thus: K(t, u) = 2e-(tu). It can easily be shown that xo(t), n(t), and K(t, u) are L2 functions, in the interval 0 < t < T, where T is any finite period of interest. Therefore, we can conclude that there exists a unique solution that belongs to an L2 space in the same interval. Moreover, an upper-bound function can be obtained by using Equation 181 and 182 as follows. Let 3(t) = Ix0(t)l + In(t)l + mca(t) (181) and 0o k A m = c ) (182) 0 Then to calculate C, use Equation 106. This gives C=/3 L /22 - 2/5 e1 - (206) Also, applying Equations 103 and 104 gives a(t) = 2(1 - e-2t)1/2 (207) and A =2 -+ T + Te 2 )(208) 42

Institute of Science and Technology The University of Michigan and T is finite. Then choosing T = 1 gives A = 0.364, C = 0.3, and m =.143. Now, by substituting into Equation 181, we can conclude that the absolute value of the system solution x(t) will always remain below an L2 function 0(t) given by -2t ttien2tb1/2 x(t) < p(t)= e2t + 1/2 (e-te 4+.2(1 -2t)1/2 (209) It would be interesting to find the exact solution for the system given by Equation 196. It should be noted that the restrictions imposed by A. A. Wolf (stated in Section 1) fit closely here, so that one can apply Wolf's technique as follows: Partitioning the differential Equation 196 at the highest derivative gives dx -2t 2 d_-e - -x-x =f(t) (210) The impulse response W(t) of the left-hand side of Equation 210 is h(t), a unit step function.16 The solution x(t) can be written in the form of a power series as follows: x(t) = K + CnQn(t) (211) n=0 where the moment Qn(t) is given by Qn(t) =Si unW(t - u)du (212) nnl 1n+1 Qn(t) - n + 1(213) Therefore, oo C x(t) =K + rn n+ (214) n=0 But when t = 0 and x = -1, the result is K = -1. Substituting Equation 214 with K = -1 into Equation 210 gives CoI n -2t c n (215 z Cntn = e2t + 1 -1 + n(215) n n+ l n+I n=0 n=O n=0'6The definition of the unit step function is 1 t>O h(t) = 1/2 t = 0 t<O 43

Institute of Science and Technology The University of Michigan Now, recall that -2t (2t)n(216) n=O Substituting Equation 216 into Equation 215 yields ~ (-2t)n n (27 c Ctn f t + E n +1 tn+1 L Cn(2)tn+2 n=O n=0 n=0 n=O from which it is possible to obtain the recurrence relation: C = (-2) + - C (n-2) (218) n n. n (n-2) (2) where C(n2), explained in Reference 47, means n-2 C C (2) z r n-r-2 (219) C(n-2): +ln —r r=0 The first few coefficients in this case are: CO = 1, C1 =1, C2 = 1/2, and C = -1/6. The following solution is obtained when these coefficients are substituted into Equation 214: 2 3 x(t) = -1 + t-T + T3_ = -..:e (220) A comparison of the actual solution x(t), given in Equation 220, in found in Table I. For comparison, x(t) and f(t) are also plotted in Figure 10. Example Two As a second example, consider the nonlinear time-varying feedback control system shown in Figure 11. This can be described by the following equation: dt + a(t) x + 0.2X2 = g(t) (221) where a(t) = time-varying parameter 2 1 t~~~~~~+ ~~t (222) 0.2 g(t) = forcing function = 2 (l+t2)2 x(O+) = 1 44

Institute of Science and Technology The University of Michigan TABLE I. SOME NUMERICAL VALUES OF Ix(t)l AND 3(t) FOR THE SYSTEM GOVERNED BY EQUATION 196 t 0 0.1 0.2 0.3 -0.4 0.5 0.6 0.7 0.8 0.9 1 Ix(t)! 1 0.91 0.82 0.75 0.685 0.62 0.563 0.512 0.465 0.62 0.37 3(t) 1 0.98 0.915 0.845 0.772 0.674 0.614 0.569 0.52 0.478 0.437.-.6x(t) X.3t 0.1 2.3.4.5.6.7.8.9 1 t FIGURE 10o A PLOT OF THE SOLUTION x(t) AND THE UPPER-BOUND FUNCTION P(t) FOR THE SYSTEM GOVERNED BY EQUATION 196 Linear g(t) + X(t) 1/S + X 2t) 1 x(t) 2x2 a (t) + + x(t) FIGURE 11. REPRESENTATION OF THE SYSTEM GOVERNED BY EQUATION 221 45

Institute of Science and Technology The University of Michigan Partitioning Equation 221 at the second term gives dx 2t 1 (223) d-t + 2 x(t) = 0.2 (223) lt 1+2 (1 t The next step is to determine the impulse response W(t, u) of the left-hand side of Equation 223. This can be obtained by solving the following differential equation: dWt, u) (t, u) = b(t - u) (224) i+t The solution of Equation 224 is 2 1 + u 2 W(t,u) 2u < t (225) =0, u>t The term x0(t), as defined before, is the solution of the following linear differential equation: dx 2t 0.2 dt + 2 x(t) = 22 (226) dt 1+t (1+t2) and is subject to the boundary condition x(t) = 1. Solving for x0(t) gives ( 1- 0.2 tan t x0(t) 2 (227) Now, since 2 1 +u 2 F[t, u, x(u)] = 0.2 1 x (u) (228) Then f[t, u x(u)] 0= 2[1 + 0.2 tan-l (u)] (229) (1 + t2)(1 + u2) Applying condition 105 of the mathematical restrictions given in Section 2.2 gives n(t) = 0.2 tan-l (t)F + 04 (tan1 t2] (230).... +t2 L(t) + — 46

Institute of Science and Technology The University of Michigan By choosing a domain D = [ Ixl < 1, 0 < u < t < T] and proceeding as in the last example, it can be shown that: 0.4(1 +22 K(t, u)= + ) (231) 2 (1 + ) Recognition of x0(t), n(t), and K(t, u) as L2 functions in the interval 0 < t < T leads to the conclusion that there exists a unique solution which belongs to the L2 space in the same interval. The solution x(t) will always remain below an L2 function 1(t), which can be calculated as follows: 2 -T.4 tan-1 04 - 12 = 04 tan (t(tan + (tan 2 dt (232) (1 + t2)2 3 2(t) = 0.16 t + (233) (1 + t2)2 A2 = 0.16T t + 2/3 t3 +t5/5 A 0.16 dt (234) (1 + t where T in finite. The choice of T = 1 gives A = 0.111, C = 0.093, and m = 0.012. Substituting x0(t), n(t), m, and a(t) in Equation 181 gives the upper-bound function 1(t): (t 2 -+ 0.4 tan (t) + 0.4 [tan 2 0008 -1 t)]3,3(t) =1 2 1+04 nt) + 0.4 [an-, + 3 [tan (t)] 1 +t (235) 0.00192 (t + 2 t3 + t5/5) 1 +t2 Therefore, (t)W > 1 2 (236) 1 +t Inspection of Equation 221 reveals that x(t) = 2 is the exact solution. Hence, x(t) < P(t) all t (0 < t < 1), as it should be. I +t Example Three Remarks: In examples one and two, the interval of interest was chosen to be finite, and there was no difficulty encountered in apply the technique presented in this section. The following example, however, is intended to illustrate to the readers how this technique can be applied to problems in the infinite interval of time. In the given domain D, one should ascertain that the functions x0(t), n(t), and K(t, u), as previously defined, belong to L2 space. Moreover, inequal47

Institute of Science and Technology The University of Michigan ities 108 should be satisfied; otherwise, the technique fails to work. It is obvious that example one works only for the case in which T, as defined before, is finite; it does not work if T is infinite. This is because as T - oo, A tends to co, as is clearly shown in Equation 208'. The problem presented in example one is again presented, but modifications have been made to permit work in the infinite interval of time. The example can now be given. Consider the system whose behavior is governed by the following nonlinear differential equation: dx 2 -2t dx + x + ax =e (t = 0, x = -1) (237) where a is a constant parameter. The parameter a must be chosen in such a way that the system response remains in a given domain D, as given by the following identity: DE [IxI < d(t)= 2e-05t0 <u < t < c] (238) Partitioning Equation 237 at the second term gives dx -2t 2 + x = e - ax (239) The impulse response of the left-hand side of the above equation is given by W(t, u) = W(t -u) = e-(t-u) (240) 17 The term x0(t), as defined before, is found in example one to be -2t x0(t) = -e (241) Since F[t, u, x(t)] = a e x (242) n(t) can be found by applying hypothesis 105 of the system restrictions: n(t) = (et - e- 4t) (243) By following the same procedures used to obtain Equation 204 of example one, the weighted nonlinear function F(t, u, x) can be shown to satisfy the following condition: F(t, u, x1) - F(t, u, x2) < 2a xmax e( x) 1-2 adt) e( )x1 - x21 (244) 1'7If Ix0(t) | is not inside d(t), the analysis fails. 48

Institute of Science and Technology The University of Michigan Thus, K(t, u) can be given by the following equation: K(t, u) = 4ae-(t-u) e 0.t = 4a-t +0.5t (245) and K(t, u) can easily be shown to belong to L2 space in the given domain because 2 t 2. -t -2t( a (t) = k (t, u)du = 16 a [e - e1] (246) 0 and A= a (t) dt= 8a (247) Also, x0(t) and n(t) as defined by Equations 241 and 243 can easily be shown to be L2(0, cc). Now, the parameter a must be chosen in such a way that the upper-bound function 3(t, a), as defined by Equation 181, remains below d(t) for all t, [0 < t < cc]. Substituting for x0(t), n(t), and ao(t) in Equation 181 gives Ix(t, a) i < M(t, a)= le-2t1 + a (et _ e -4t) + m(e-t -2t)1/2 (248) where c k m = c ). (249) k=1 and c is given by 2 00 2 a c n2(t)dt = (250) 40 The parameter a, as it appears in Equation 248, controls both the absolute value of x(t) and the upper-bound function 1(t) in such a way that inequalities 108 must be satisfied. This, of course, permits flexibility in the application of the given technique. If the parameter a were not present, the inequalities 108 would not be satisfied and the technique would fail to work.18 8Where this is true, 1(t) may assume the shape of the curve /*(t), shown in Figure 12. 49

Institute of Science and Technology The University of Michigan An investigation of Equation 248 reveals that a value of 0.25 assigned to the parameter a satisfies inequalities 108 and gives the following results: ca(t) (e-t _ e2t)1/i2 n(t) = 12 (e-t e-4t) c = 0.039 A = 0.707 m = 0.048 Substituting from Equation 251 into Equation 248 gives Ix(t)I <:(t) = e- 2t + 0.083 (et - e 4t) + 0.048 (et _ e-2t)1/2 (252) The upper-bound function P(t) and the function d(t) are plotted in Figure 12 for the purpose of comparison. Some numerical values of these two functions can also be found in Table II. I T 2 1.8 d(t) 1.6 1.4 1.2 - *(t) 1\ 1' ( \\\ 0.8 \ 0.6 /3(t) C 0.4 0.2 1 2 3 4 5 FIGURE 12. A PLOT OF d(t) AND P(t) FOR THE SYSTEM GOVERNED BY EQUATION 237 50

Institute ot Science and Technology The University of Michigan TABLE II. SOME NUMERICAL VALUES OF d(t) AND i(t) FOR THE SYSTEM DESCRIBED IN EQUATION 237 t 0 1 2 3 4 5 co d(t) 2 1.216 0.7220 0.4464 0.2706 0.1638 0 O(t) 1 0.1989 0.0447 0.0171 0.00832 0.00451 0 To obtain the approximate solution of the new system equation, namely: dx x -2t dt rx +-4 = e (253) one has to substitute for C and A in the error Equation 170 2 2n cc k, 2 C A A (170) k=0 Then, continuing the above substitution in Equation 170 gives 6 2 (0.038)2 (2.5522)2 (0.5)2 (254) n (n - 1) (254) For n = 1, Equation 254 gives 6 2 = 0.0039 (255) The first-order approximate solution x1(t) for the above value of mean square error 61 can be given by xl(t) = x0(t) - n(t) (256) Substituting from Equations 241 and 243 into Equation 256 gives x1(t)= -[e-2t - (e - e4t)] (257) Inspection of Equations 257 and 252 reveals that the modulus of the approximate solution xl(t) remains below the upper-bound function 3(t) for all t, as it should be. 51

Institute ot Science and Technology The University of Michigan 2.10. CONCLUSIONS AND REMARKS The partitioning technique has been used to analyze a certain class of nonlinear systems whose dynamic behavior can be represented by a nonlinear differential equation. It was found that when certain restrictions are placed on the linear, nonlinear, and forcing function terms, this equation has a unique solution which is the limit of a sequence of iterates {xn}, each of which, together with the system output x(t), belongs to an L2 space. The solution was shown to remain below some upper-bound function:(t) in the L2 space for every t in the period of interest 0 < t < T. The upper-bound function:(t) defined by Equation 181 can be easily calculated since it is a linear combination of functions defined by the class of system under consideration, and it may be changed by changing the system parameters. This method should help control engineers gain insight into the stability and boundedness of system responses, and hence, is a first step toward solving the problem of synthesis. Equation 170, derived in this section, gives an upper bound for the error when the actual solution x(t) is approximated by the n-th iterate xn in the sense of the mean square. The number of iterations n required for a given approximation was shown to depend on how the system equation is partitioned and, therefore, one can compare, with the help of Equation 170, different partitions in order to have n as small as possible. It should be pointed out that since the L2 space isacomplete space, the solution x(t) is guaranteed to be in the same L2 space of the sequence of iterates {xn}, a fact which is not true for the case studied by A. A. Wolf in Reference 47. This is because the limiting process of a sequence of functions of exponential type may give rise to functions of higher order but not of the same exponential type. The following example, introduced by A. A. Wolf [47, p. 22], illustrates this: 2 -t x +2 et (258) xo(t) = - et (259) x (t) (1 - et) i x (u)du (260) then, as n - oc, xn approaches infinity. The resulting limit is clearly of higher order and not of exponential type. In a given domain D, the modified Lipschitz condition defined by Equation 102 is sufficient, but not necessary, for the existence of Picard iterates txn} and the uniqueness of the solution x(t). The Lipschitz condition implies continuity of the weighted nonlinear function F[t, u, x(t)] with respect to x but not to t or u. The physical significance of the condition is that when t and u are 52

Institute of Science and Technology The University of Michigan fixed, the resulting function F(x) has bounded difference quotients; the bound is given by the L2 function K(t, u) defined in this section [29, p. 38]. If the domain D is convex (that is, D contains the line segments connecting any two points in D), then an application of the mean-value' theorem of differential calculus will show that the existence and boundedness of the partial derivatives of F(t, u, x), with respect to (x) in the given domain, are sufficient for F(t, u, x) to satisfy the Lipschitz condition [8, p. 54]. This means that the rate of change of the weighted nonlinear function should not be rapid. Nonlinearities with sharp breaks would not satisfy the Lipschitz condition and should be avoided, but since they are not too likely to occur in nature, the given theory can be applied to a broad class of physical systems [47, p. 191]. Although uniqueness does not always imply the convergence of successive approximations, the restrictions placed on the theorems given in this chapter are sufficient to assure convergence [8]. (Examples of convergence are given in Appendix E.) The conversion of the system differential Equation 96 into the equivalent integral Equation 122 proves an interesting relation between circuit configuration and feedback-system configuration. Equation 96 is composed of both a linear and nonlinear set of terms and a forcing function; therefore, it can represent the electric circuit of Figure 13. The integral Equation 122 can represent the feedback system of Figure 14, which has the same output as the current flowing in the circuit of Figure 13. This establishes a relation between circuits and feedback systems. The class of forcing functions g(t) allowed in this work is defined completely by condition 102. This condition states that the forcing function should be such that when applied to a linear system consisting.of the linear terms L(D, t) (defined by Equation 97 and shown in Figure 15 the energy associated with the response x0 of the linear system should be finite.19 This condition includes many physical systems of interest. By an application of the Schwarz inequality to Equation 99, it can be shown that llx011 < iigll 11w t, ull (261) i.e., condition 102 can be satisfied by choosing functions g(t) and W(t, u) that have bounded norms. Finally, it is worthwhile to note that the method of analysis presented in this chapter is best suited for handling nonlinear systems with parameters, because it gives control engineers an easy and flexible means for choosing parameters that will satisfy certain requirements in operation during a specified time. This is clearly shown in example three and will be elaborated upon in the following section.'9See Reference 25, page 123, for the following definition: "For an arbitrary real valued function of time x(t), we call the quantity x(t)2 the instantaneous power associated with x(t). The total energy associated with this function is the integral J x(t) dt when this integral exists." 53

Institute of Science and Technology The University of Michigan -[ L(D, t) Linear t) Element N Nonlinear Element FIGURE 13. ELECTRIC CIRCUIT CONFIGURATION OF EQUATION 96 Input g(t) +x(t) = Output L(D, t) N[x(t)] (t) N(x) x FIGURE 14. A FEEDBACK SYSTEM REPRESENTATION OF THE NONLINEAR INTEGRAL EQUATION 122 g(t) L(D, t) x0(t) FIGURE 15. A LINEAR SYSTEM REPRESENTATION OF EQUATION 97 54

Institute ot Science and Technology The University of Michigan 3 STABILITY ANALYSIS This section discusses the use of the partitioning technique and the state-variable technique to study a class of control systems whose behavior can be represented by Equation 96. Placing suitable restrictions on the system equation revealed that the state variables (v1, v2,..., vn), which represent the state of the system in the phase space, belong to an L2 space. These restrictions include those presented in Section 2 as sufficient, but not necessary, for the existence and uniqueness of the system's response. Two definitions of the norm of the state of the system were used to study asymptotic stability in the Lyapunov sense. This method is useful for an analysis of the system trajectory at any time t in the period of interest. An expression for the upper-bound state is calculated, and an example is included to illustrate the method presented. 3.1. SYSTEM RESTRICTIONS The notations, definitions, and system description introduced in Section 2 are followed in the present analysis. The system must satisfy a number of mathematical restrictions which define the class of admissible systems under consideration. These restrictions are: (1) x0(m)(t) belong to an L2 space over an interval [0, T], and are assumed bounded. m-= 0, 1,..., (n - 1), where n is the order of the highest derivative of the linear part of the system equation. (m) = the m-th order derivative with respect to the independent variable t. (262) (2) In a given region D, [I xI< d(t), 0 < u < t < T], F(m)(t, u, x)20 satisfies the following two conditions: (a) F(m)[t, u, x(u)] satisfies the modified Lipschitz condition; that is, if the two triples (t, u, Zl) and (t, u, z2) are in D, then I F(m)(t, u, Z) F(m)(t, u, z2)[ Km(t, u,)[z I - 2(263) where [Km(t, u)] are L2[O, T], that is K (t, u) du < a (t) T (264) a m(t) dt < A 20 ) 1~ amF(t, u, x) F(m)(t, u, x) F(tu, x) and, as before, F is a function of x only, and not x, etc. atm 55

Institute of Science and Technology The University of Michigan (b) F(m)[t, u, x0(u)] du < n m(t) where nm(t) are bounded and belong to L2(0, T), that is, (t) 2 where () and m are as defined above(266) where (m) and m are as defined above 2 2 A and C are finite positive numbers m m T is the period of interest (0 < T < + oo) The functions a 2(t), nm(t) are continuous and bounded over [0, T] (3) Condition 108 is still required; that is, q(t) < d(t) 1 x0(t) I < d(t) (267) where x0(t) and q(t) are as defined by Equations 99 and 107. (4) The time-varying parameters a,(t) and I N(x, t) - g(t) are bounded over [0, T]. (268) It should be noted that when m = 0, restrictions 262 through 268 include the set of restrictions defined in Section 2.2. These restrictions, that is, those in which m = 0, were shown to be sufficient for the existence and uniqueness of the solution x(t), which is the limit of a sequence of iterates {x(t)} belonging to the L2 space, and given by the following equation: t x+1 = xo(t- F[t, u, xn(u)] du (n = 0, 1, 2,...) (269) 3.2. THE STATE-VARIABLE REPRESENTATION Consider the system described by Equation 96 under the restrictions placed on it by Equations 262 through 265 and let: x =V1 dx (270) dt V n-1 dtn-l n 56

Institute of Science and Technology The University of Michigan These equations become dv dt =v2 =fl(t, V1v' v2'' Vn) dv dt2 3 =f2(t, v1, v2,, vn) (271) dv dt an() 1 an-l(t) v2'* + N1 =fn(t, Vv12 n) where N1 = N1(x, t) = g - N (272) Equation 96 can be put into the following matrix equation form: v O 1 0 v V2 0 0 0 v2 0 = + (273)'t 0 0.. Vn_ a-a(t)-a_ (t) -a(t) V N Again, Equation 273 can be represented by the following two equations: d - d [v] = [A(t)] [v] + [N1] (274) d [v] =f(t, v) (275) where v = a vector comprised of the state variable (v1, v2,..., vn); it can be used to represent the state of the system in the phase space in which each succeeding axis represents the rate of change of the quantity measured along the one preceding it f = a column vector made up of the functions (f1, f2''...' fn) as defined by Equation 271 [A(t)] = an n by n matrix, as defined by Equation 273 57

Institute of Science and Technology The Uhniversity of Michigan 3.3. THE PROPERTIES OF THE STATE VARIABLES (v1, v2,..., vn) OF THE CLASS OF SYSTEMS UNDER CONSIDERATION (THEOREM E) Theorem E summarizes some of the properties of the state variables which are of interest in this study. Theorem E Under the restrictions (given by Equations 262 through 268) placed on the linear, nonlinear, and forcing function terms, the state variables vi(t) belong to L2[0, T], and are bounded by corresponding arbitrary L2[0, T] functions 0i(t). Since INl(t, vl)l and lai(t) I, which were defined before, are assumed bounded over the interval [0, T], then the state variables vi(t) are uniformly continuous21 over the same interval [0, T]. In particular, if T = +oo, then it follows that Lim vi(t) =0 (i = 1, 2,...,n). Proof The following two lemmas are needed to prove the theorem: Lemma One22 If a function f(t) has a bounded derivative on any interval (finite or infinite), then f(t) is uniformly continuous over the same interval. 28 Lemma Two [50, p. 86] If a function f(t) is uniformly continuous over the interval [0, oo] and if f2(t) dt is finite, then it follows that Lim f(t) = 0. t-oc Recall from Equation 122 that the state variable v1(t) = x(t) can be given by x(t) = v1(t) = x(t) - F[t, u, v1(u)]du (276) and that t t dt g(t, u) du = g(t, t) + g(t,u)dt (277) J O J 0 21See Appendix H for the definitions of uniform continuity. 22Lemma One can be easily obtained by an application of the law of the mean: f(t2) - f(t) I = f(t*) lt2 - tl < KIt2 -t11, t1 t* t2 23 See also Appendix H. 58

Institute of Science and Technology The University of Michigan Upon differentiating Equation 276 and using Equation 277, the following equation is obtained: dx dv1 (t) d(278) dt - dt dt x0(t) + F[t, tvl(t)] - dt F[t, u, Vl(U)]du However, since W(t, t) = 0, it follows that F[t, tvl(t)] = W(t, t) N[V1(t)] = 0 (279) Equation 278 can be put in the form: x(t)() =v2 =v1(1) = x()(t) - ()t, u, v (u)]du (280) where, as defined before, the superscript (1) denotes the differential operator d/dt, and F( )[t, u, v1(u)] = W(1)(t, u) N1[vl(u)] (281) From the properties of the impulse response,24 it follows that W(t, t) = W( )(t, t) =.. = W(-2 )(t, t) = 0 (282) W(n 1)(t, t) = 1 a (t) Using Equation 282 and differentiating Equation 276 k times yields x(k)= Vk+1 = x(k)(t) - | F(k)[t, u, v1(u)]du (283) Equation 283 gives all the state variables, v1, v2,..., vn. When k = 0, it yields the state variable v1, which was shown by theorems A, B, and C to exist, to be unique, and to belong to an L2 space. Now, to prove that all the vi, (i = 1, 2,..., n) in general belong to an L2 space, it is necessary to define a sequence xm according to the following recurrence relation: (k) x0(k)(t (k) m+1) 0 )(t) - F(k)[t, u, x (u)]du (284) where k = 0, 1,2,...,.n - 1 m = 0, 1, 2,... 24See Appendix C. 59

Institute of Science and Technology The Uhniversity of Michigan Equation 284 can be used as follows to obtain xm(k): m xm() = ( (t) - F( )[t, u x 1(u)]du (285) Subtracting Equation 285 from Equation 284 and squaring gives x (k) (k))2 ( {F(k)[t u x (u)] - F(k)[t u x (u)] du) < Kk(t U) Ixm - Xm-l < Kk(t u)2 du [Xm(u) - x m(u)]2du (286) Equation 264 can be used to put Equation 286 into the following form:26 (X1 (k) (k))2 < ak(t) [Xm(u) - xml(u)] du (287) 0 But Equation 134 gives 2 m-2 T C A x [m () -m 1 (u)]2 du < (m (288) Therefore, substituting Equation 288 into Equation 287 gives 2 C2A m-2 +(k) - x (k 0 0 2(t) (289) m+1 mi (m - 2) k where m = 2, 3,... k = 0, 1,.., n - 1, n = the order of the system Equation 289 shows that all m(k)m are L2 functions. It can easily be proved that Lim x ( exists and is given by26 (k) (k) Lim x = x v (290) m-co m k+1 where k = 0, 1,..., n - 1, n= order of the system. 2C [Xm(U)\- xm1(u)]2du< 0 d[Xm(U-, Xm-l()]2 d' t < T 26See Appendix G. 60

Institute of Science and Technology The University of Michigan Therefore, it can be concluded that all the state variables vi (i = 1, 2,.., n) belong to an L2 space in the assumed interval of interest T (0 < t < T). Using Equation 528 of Appendix G, an upper-bound function for each of the state variables vi (t) can be calculated so that vi(t) < i(t) for (i = 1, 2,...,n) (0 < t < T) (291) The functions.i(t) satisfying inequality 291 are given by ~i(t) = x0(i-1)(t) + Ini 1(t) i + m oi (t) (292) where x0( (t), ni_1(t), and aci (t) are functions defined by the restrictions placed on the system for i = 1, 2,..., n, and m is a positive constant given by m = Cx AAk (293) k=0 It is clear from Equation 292 that the functions pi(t) lie in the L2 space, since the right-hand side of Equation 291 is a linear combination of L2 functions. Also, the upper-bound functions can be changed as desired since they depend only on known functions defined by the system restrictions. If the period of interest T is chosen to be infinite, then from the properties of the L2 spaces, we can conclude that v.2 dt < o (i=1,2,...,n) (294) It is not necessary for a.continuous function which is square integrable on the interval [0, +cc) to have a zero limit at t = +Co since in general the limit may not exist. However, by considering the matrix Equation 273, it can easily be shown (since the time-varying parameters a (t) n and the nonlinear function N (t) are assumed bounded) that each of the state variables vi(t) has a bounded derivative in the interval of interest (assumed infinite). Thus Lemma One shows that vi(t) for i = 1, 2,..., n are uniformly continuous over the same interval. Now, by considering Equation 294 and by applying Lemma Two, it can be concluded that Lim v.(t) = 0 for i = 1, 2,..,n (295) t-oo 1 27See Appendix H. 61

Institute of Science and Technology The University of Michigan 3.4. STABILITY CONSIDERATIONS Consider the system described by the vector Equation 275. An equilibrium position ve exists such that all the state derivatives dvi/dt for i = 1, 2,..., n are simultaneously zero. It can easily be seen from Equations 271 that, for the system under consideration, the origin of the phase space v = (v1, v2,..., vn) = 0 is a position of equilibrium if the condition N1(V, t) = 0 at v(t) = 0 (296) is satisfied. Thus, under condition 296 the solution vr = 0 is an equilibrium position of Equation 275 and its stability is in question. With the passage of time, the state vector v whose components are v, v2,..., vn traces a curve in the phase space known as the system trajectory which describes all possible states of the system. Figure 16 is a possible trajectory in the n-dimensional space. It may be thoughtof as the projection of a solution of the equation, plotted in a space containing the n-space axis and time axis, into the state space. This is shown in Figure 17. v iV /'v~~~V2 V1 FIGURE 16. SYSTEM TRAJECTORY 62

Institute of Science and Technology The University of Michigan V2 Traj ectory Solution v1 t FIGURE 17. PROJECTION OF A SOLUTION CURVE [12] A simple measure of the departure of one possible state i from the equilibrium-position v = O can be obtained from the norm of v denoted by v1 i1. The norm is defined as a function which assigns to every vector v in a linear space a real number denoted by KIvI so that [1, p. 406]: vI l is defined for every v in the space ivll = 0 if and only if v = O lvil > Ofor all v > O (297) Ilu+ vt < Ulvll + lu-li for allu and v a1v l = Ial IV -i1v for all v, and a is a real constant Some commonly used norms are [8, p. 17]: V.; (' i5(298) lvi = max {v.} (299) n lvll = Y Iv. (300) i=l1 63

Institute of Science and Technology The University of Michigan 1/2 1111 f v(,.., V(Sl s ds, dsj (301) where R is a fixed region of the s1,..., s space, and Ivi denotes the modulus of the vector v in the usual way. It can be shown that the norms, as defined by Equations 298 through 301, satisfy the conditions for the norm as stated under definition 297. Thus, it can be said that for the system under consideration that an output state u is greater than an output state v if and only if If a system is perturbed slightly from its equilibrium state at the origin, it is desirable in many applications to have all subsequent motions remain in a correspondingly small neighborhood about the origin. This is shown diagramatically in Figure 18 and means that the system under consideration is stable. However, in order to investigate stability it is necessary to define precisely the concept of stability to be applied. V2 FIUC v(t) FIGURE 18. GEOMETRIC DEFINITION OF STABILITY. C = upper-bound state, V = system trajectory, V0 = initial state, O = equilibrium state. Definition of Stability [7, 12, 14] The equilibrium position v = 0 is called stable in the sense of Lyapunov if, for every e > 0, there exists a 6(E) > 0 so that vIv(t)I <I e whenever I V(0)1 < 6 for all t > 0. In other words, the equilibrium position v = 0 is stable if the initial magnitude of the norm of the state vector is small enough so that 11 vii can be made to remain below an arbitrary upper bound. If v(t) approaches zero as t approaches infinity, that is, if Lim v.(t) = 0 (i = 1, 2,.., n) t-o 1 64

Institute of Science and Technology The University of Michigan the equilibrium position is called asymptotically stable. If 6 is independent of the choice of the initial time, then it follows that the system is uniformly asymptotically stable. It is necessary to confine the application of this definition of stability to those cases which satisfy the definitions of the norm ivI given by Equations 297. However, for the purposes of this report, there are two reasons which make it desirable to study the behavior of the system under the two definitions of the norm given by Equations 298 and 301: (1) The Euclidean norm defined by Equation 298 gives the magnitude or length of the state vector v(t) during the period of interest. (2) The norm defined by Equation 301 gives a measure of the mean square value of the state vector v(t) during the period of interest. Hence, the two definitions of the norm lv I adopted in this study govern the transient behavior of the system under consideration, at any particular time t and also on the average during a given interval of time, a situation which is of interest to many system engineers. Now, if the functions x0(m)(t), for any initial state v0, satisfy conditions 262 of the restric(m) tions placed on the system, i.e., if x0 (t) belongs to an L2 space, and if conditions 263 through 268 are true, then it can be concluded that (1) The system under consideration possesses a unique solution existing in a domain to the right of the initial state. (2) The system trajectory v(t) approaches 0 as the time t approaches infinity. (3) The norm IvlI, as defined by Equation 298, is bounded by an arbitrary upper norm l:(t) l given by tI = l i (t 11j~t~ll ii kt~~''" (302) where hi(t) is as previously defined by Equation 292. (4) The norm v Iv as defined by Equation 303 (a special case28 of Equation301), namely [IvI =j I v(t) 2 dt (303) exists and is less than an arbitrary positive number N, given by the following equation: N2=J (t)c 2dt (304) For the special case, v = 1, s = t, and R = the interval 0 < t < oc. 65

Institute of Science and Technology The University of Michigan where [-(t)l = x2(tj/ (305) The existence of the two integrals given by Equations 303 and 304 can easily be seen because both of the vector functions v(t) = [vl(t),..., vn(t)] and.(t) = [.l(t),..., n(t)] are L2 functions, as previously proved. (5) The system trajectory v is uniformly continuous with respect to the initial conditions.29 Consequently, the system under consideration is asymptotically stable according to the definition of stability given previously. For autonomous systems which are free and stationary, asymptotic stability implies uniform asymptotic stability. This is not true for nonstationary systems which are, in general, the systems discussed in this report. Unfortunately, in nonstationary systems, asymptotic stability alone does not guarantee that bounded inputs will give rise to bounded outputs; uniform asymptotic stability is needed. This is indicated by the following example suggested by Kalman and Betram [14]: Consider the nonstationary system dx + = g(t) 0 > t (306> t) (For values of to less than 0, this equation has finite escape time and does not define a dynamic system.) The impulse response of this system is given by to W(t, t0) = t t0 > 0 (307) Evidently, the system is asymptotically stable. However, the unit step response of this system, which corresponds to a zero state at some time t > 0, is given by x(t) = ) dT = (t, - t)d =/2t (308) which tends to co with t' 29This statement is deduced from arguments presented in Section 2.8. 66

Institute of Science and Technology The University of Michigan 3.5. AN EXAMPLE CALCULATION: THE ANALYSIS OF A SECOND-ORDER NONLINEAR CONTROL SYSTEM Consider the system equations: dv dt = 0(V +12v2' ) (309) dv2 dt = K21v + K22v2 + K2N(t)N(vl' v2) = f2(t' 2' t) with initial conditions at t = 0, v1 = v2 = 0 where K2N(t) = a time-varying parameter K12, K22v K21 = constant parameters v1 and v2 = the state variables that represent the state of the system in the phase plane at any instant of time N(v1, v2) = a nonlinear term Combining Equations 309 by eliminating v2 gives 2 d v1 dv dv -- K v+K 1 1'K t) Nv (310) K12 dt 21l + K22 dt + K2Nt) l, K1dt Equation 310 is equivalent to the following equation: v1 = K21K12V1 + K22Vl + K12K2N(t) N V K (311)2 12 dv1. d v where v - vandv = 1 dt 1 2 dt Consider the nonlinear control system shown in Figure 19. The system shown in Figure 19 can be described as follows: Jc(t) + (B + KaKb)t) + KaKcc(t) = KaK2N(t)N(v) (312) Recall that v1 = r(t) - c(t) 1 = r(t) - c(t) (313) v1 = r(t) - c(t) 67

Institute of Science and Technology The Undiversity of Michigan K2N(t) r (t).XX K C(t) N(Vl, v2)a Js + Bs KC FIGURE 19. A SECOND-ORDER NONLINEAR CONTROL SYSTEM Equation 312 gives J1!l(t) + (B + KaKb)Vl(t) + KaKcVl(t) + KaK2N(t) N(Vl) = Jr(t) + (B + K Kb)t + K Kcr (314) To study the homogeneous differential equation of the system Equation 314, let r(t) = O. This gives B +KK K K K vl(t + J +v) + ( t) + ( )K2N(t) N(1)= 0 (315) The homogeneous Equation 315 of the second-order control system shown in Figure 19 is identical to Equations 309 or their equivalent Equation 311 if B +K K a b = -K22 KK a c=-K21K12 (316) K a -= -K12 68

Institute of Science and Technology The University of Michigan Equations 316 relate the parameters K12, K21, and K22 of Equations 309 to the parameters B, Ka, Kb, Kc, and J of the control system. This means that changing the set of parameters in Equations 309 is equivalent to changing the control system parameters; the converse is'also true. From Equations 309 it can be seen that if N(vl, v2) = 0 (317) at v1 = v2 = 0, then the origin of the phase space will be an equilibrium position at which dvl/dt and dv2/dt will be simultaneously zero for all t. Let the functions of interest for N(vl, v2) and K2N(t), be N(v1, v2)= V1 (318) +2at K2N(t) = e where a is any arbitrary constant, yet to be determined. Now, the problem at hand is: Given a nonlinear function N(vl, v2) and the specific form of K2N(t) for the control system shown in Figure 19, choose the parameters K12, K21, K22 and a so that (1) The actuating error signal v1 and the rate of change of signal vl belong to an L2 space in the interval 0 < t < c. (2) Given two specific numbers H1 and H2, the system satisfies the following requirements: Hi 2' 0 vl (t)dt < H1 (319) 1 (idt < H2 This means that the unknown parameters are required to keep the integrated square error [vl(t)] and the integrated square rate of change of the error [vl(t)] below certain specified levels. Equations 319 also require that the state variables vl and v2 = vl should finally reach the equilibrium position at the origin of the phase plane, shown in Figure 20. To solve the problem the required parameters should be chosen so that the system under consideration satisfies the restrictions already studied. Now, x0(t) is the solution of the linear equation: B+ KbK KK vl(t) + J Vl(t) = o (320) with the initial conditions that t = O, x0 = 1, aId x = 0. 69

Institute of Science and Technology The Oniversity of Michigan V1 I.~~~~ t v (a) (1, 0) V7~so (b) FIGURE 20. REQUIRED SYSTEM CHARACTERISTICS. (a) Time domain. (b)Phase plane domain. Therefore x0(t) and x0(i)(t) take the form 1 mIt(1) x0(t) - 2 e + e X (t) =m m2 -mi m2 emi + m1m2 e) where mi and m2 are the roots of the characteristic Equation 320 as follows: B + KaKb (B + KKb) K Kc mi = - 2J. (i = 1 and 2) (322) 2J + 2 J Also, W(t, u) is the solution of the following equation: B+KK K K vi(t) + J vl(t) + J vl(t) = 6(t, u) (323) 70

Institute of Science and Technology The University of Michigan with the initial conditions that at t = u, W(t, t) = O, W( )(t, t) = 1 where 6(t, u) is the unit impulse applied at t - u. Therefore, rml(t -u) m2(t-u) 1 Wt ) W(t - u ) m1 -[e m 1 e m(tu h(t - u) (324) where mi1, m2 are as defined before, and h(t) is the unit step function. The weighted nonlinear function F[t, u, x(t)] can now be given: F(t, u, v1) = mm3 [m m1(t-u) m2(t-u)] +2au2 (325) F (t, u, v e - e v1 (325) where m3= K /J (326) Also, (1 m F ml(t-u) m(t-u +2au 2 F( (t,u, v1)= m1 e - e2 e+2au 2 (327) 1-m 1 (3272 Both F(t, u, v1) and F(1)(t, u, v1) satisfy the Lipschitz condition in the domain D= [|x < 1 (O < u < t < oo)], since 2m3 [emlt (ml-2a)u m2t -(m2-2a)u] IF(t, u, x1) - F(t, u, x2) < m- e - e e x 1 2 (328) 2me 2m3 [ ml(t) -(m-2a)u m2(tt)-(m2-2a)u]l F()(t, u X1) - F( (t, u, 2)1 < m) - m 1 2e - e e x21(329) Therefore, from Equations 328 and 329, we have 2m3 [emlt -(m+.2a)u m2t (m2-2a)u K0(t, u) m1 - m2 e- e e -( 2-2a)u] 2m3 [ mt -(m1-2a)u m2t -(m2- 2a)u] K1(t, u) = m [ m2 e1 e - 2 e (331) 1 m2 22e (3371 71

Institute of Science and Technology The University of Michigan Also, applying Equations 265 of the restrictions gives t m3 ml(t-u) m2(t- u) m lumu n (t) = { e - e -m2 e + m1 e e du (m1 - ) 2(ml+a)t 2(m2+a)t (m1+m2+2a)t mit m t = C01 e + C02 e + C03 e + C04 e + C e (332) and 3 ml (t-u) m2(t-u) u mu2 +2audu n1(t) = efA [3 m - m2 e m2 e + ml e e (mI - m2) 2(ml+a)t 2(m2+a)t (ml+m2+a)t mit m2t Cll e + C12 e + C13 e + C e + C15 e (333) w-+ 2 2m2 + 2a - m2 m3. 2 _ 2 where C i _ _ i 2 -1 ( 2in 2 m 2 1m =2 22) in3 1in2 2i m3 ml ml,3I m1 - m2 (m- m2) +2a 2 m2 + m3' I-m2 2mm2 m12m C04 _ m + 2 + 2 m2 + 2- a1 (mi1 - m2 2) m3 m 2 2mm m2 m C05-(m ( m 3 m + 2 22a 2 2+ 2a - ml -05 (m1 m)2mL +2- mi2 m1 + 2a + im2 + 2a 72

Institute of Science and Technology The University of Michigan (Equation 334 continued)3 2 m3 ___ _ 1 1 m C12 m 3 2m2 + 2a - r m2 + 2a (m- m2r2) 2 2 m3 -2m 22m2mm2 C 3_ _ 1 2 1 2 1 (m3 = m2 + +2 (mi+ 2a)J C 1 2m C14 = 1mC04 C15= m2C05 Applying condition 264 of the restrictions gives 2 t2 (m1+m2)t 2m1t 2m2t d +4at a0 (t) = Ko (t, ) du =d01e + d02 e +d03 + d e d4 (335) and aI;otK2 = 2fl (m1 +m2)t 2m t 2m2 t a= K1 (t, u)du=d e +d e d1 e + d14 e (336) where 4m 2 01 = 2 mI m2 - 4c (m1 - m2) 1 2 4m 2 3 1 4m 2 3 1 03 (m - m2)2 2(m2 - 2a) 4m3 1 ( 1 2 d + (337) d04 - 2) +2(m - 2 2 - 2) (ml + m2 4a) (in1 in2) 1(ml dll 1 mlm2d01 d12 = m12d02 73

Institute ot Science and Technology The University of Mich ig a n 13 2 d03 2 2 2 -4m3 m1 m2 2 14 (m - m2)2 2(m1 - 2a) 2(m2 -2a) - m + m2 - 2a An examination of the functions x 0(i)(t), ni(t), and ai(t) for i = 0 and 1 reveals that these functions belong to an L2 space (in 0 < t S< o), if the values of ml, m2, and a (hence, the system parameters) are properly chosen. These functions are also uniformly continuous and, if they belong to L2 space in the infinite interval, they will have zero limits at t = + c. In order to evaluate the two integrals in Equations 319, one has to make use of Equations 292 and 293 to evaluate A3 and CO, as follows:30 A0 = LI 02(t) dt 2 2m 4a (338) 00 1/2 - 2 2 C2 2 271/2 C 2(t) dt] = 4 02 03 04 05 (339) C~4(m o+ a) 4(m2 a) -2(ml +m 2a) m Cm C0 n~~(t) = - +C2(m satisfies the requirements: a = Kb = 1 B = K 2 (340) a= 1.5 The parameters of Equation 340 correspond to m1 = -1, m2 = -2, m3 = 1, and a = -1.5. 1 2 K22 =-3 (341) K 1 = +2 30 In evaluating these integrals, it was assumed that mi1, m2, and a are negative numbers. Use was made of the fact that (a + b + c + d + e)2 < 4(a2 + b2 + c2 + 2d2 + 2e2). 74

Institute of Science and Technology The University of Michigan By using the set of parameters derived in the previous analysis, the Runge-Kutta method, and a digital computer, the system was analyzed numerically. The values of v1(t) and v2(t) for different time intervals are shown in Table III and are plotted in Figure 21. The phase. plane behavior of the system is shown in Figure 22. Notice that this type of analysis in which the interval of interest is infinite is not suitable for digital computation from a practical point of TABLE III. SOME NUMERICAL VALUES OF v1(t), v2(t) FOR THE SECOND-ORDER NONLINEAR CONTROL SYSTEM t 0.1 0.2 0.3 0.4 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 v1(t) 0.99 0.954 0.908 0.857 0.802 0.538 0.344 0.215 0.134 0.081 0.049 0.027 0.017 0.010 v2(t) -0.245 -0.401 -0.493 -0.541 -0.556 -0.464 -0.318 -0.205 -0.129 -0.080 -0.049 -0.027 -0.016 -0.010 1.0.50 0 0.5 1.0 1.5 2.0 2.5 3.0 - __" " ___ -.25 -.50 -.75 v2(t) FIGURE 21. THE SYSTEM TRANSIENT RESPONSE 75

Institute of Science and Technology The University of Michigan V2 0.1.2.3.4.5.6.7.8.9 1 Vi -.1 -.2 -.3 -.4 -.5 -.6 - FIGURE 22. PHASE PLANE TRAJECTORY OF THE SECOND-ORDER NONLINEAR CONTROL SYSTEM view. However, in the example chosen the values of v1(t) and v2(t) happened to decay very quickly with time so that beyond a certain time interval the response functions v1(t) and v2(t) were considered approximately equal to zero. 3.6. CONCLUSIONS AND REMARKS The partitioning technique and the state variable approach have been used to analyze and study a system whose dynamic behavior can be represented by a nonlinear differential equation containing some linear terms, some nonlinear terms, and a forcing function term. Section 2 proved that this system possesses a unique solution which belongs to an L2 space, and Section 3 proved that under suitable restrictions all the state variables vl, v2,..., vn belong to an L2 space. By studying the definitions of asymptotic stability in the sense of Lyapunov, the system, as restricted, was found to satisfy the definition. Hence, the restrictions placed on the system are sufficient, not only to prove uniqueness of a solution which belongs to L2 space, but also to produce asymptotically stable systems in the sense of Lyapunov. Further, by selecting two different definitions of the norm, it was possible to study the behavior of the system trajectory at any time t and on the average during the interval of interest. Studying the average behavior of the system trajectory is not usually desirable because in many applications, although the average value of the system trajectory cannot exceed a certain specified value, the actual sys76

Institute of Science and Technology The University of Michigan tem response may have an undesirable shape as shown in Figure 23. Therefore, an expression is obtained for an upper-bound state vector within which the system always remains during operation. This upper-bound state vector is found to depend only on the system restrictions; hence, it can be changed to suit a specific application. An example was given to illustrate the method presented. The author believes that the analysis presented herein should help system engineers to analyze or synthesize systems containing nonlinear characteristics. v(t) If r~~Average Value 0 t t=T FIGURE 23. ACTUAL AND AVERAGE SYSTEM BEHAVIOR 4 EXTENSION OF THE ANALYSIS TO MULTILOOP SYSTEMS This section extends the analysis to multiloop nonlinear systems whose dynamic performance can be represented by n simultaneous nonlinear differential equations with time-varying coefficients. The procedures are similar to those used in the analysis presented in Section 3. The new admissible class of systems is defined by placing suitable restrictions on the linear, nonlinear, and forcing function terms. These restrictions are found to be sufficient, but not necessary, to produce unique and stable multiple outputs all of which belong to an L2 space in the interval of interest. The multiple outputs are given as the limit in the mean and, in the ordinary sense, of a sequence of functions in the L2 space. 4.1. INTRODUCTION Thus far we have analyzed a class of systems whose behavior is described by a single equation (96). That equation can be considered the canonical form for nonlinear systems since any 77

Institute of Science and Technology The University of Michigan such system can be described by an equation of that type. However, there are situations in which: (1) The system equation will not reduce easily to the canonic form. (2) A system contains coupling between different subsystems. (3) A system has more than one input (a multiple-input system). The input may be desired or undesired; for example, noise within the system would be undesired. The noise should be a known function of time. (4) A system has more than one desirable output (a multiple-output system). 4.2. SYSTEM DESCRIPTION AND DEFINITIONS It is of interest to control engineers to extend the analysis to include the situations cited. Therefore, it is necessary to consider the general case of n simultaneous nonlinear differential equations of the following form: L11(D, t)xl(t) + N11[Xl(t)] + L12(D, t)x2(t) +... + Lln(D, t)xn(t)+ Nln[x(t)] = gl(t) L 21(D, t)xl(t) + N21[x1(t)] + L22(D, t)x2(t) +... + L2n(D, t)xn(t) + N2n[x(t)] = g2(t) (342) L n(D, t)xl(t)+ Nn2[xl(t)] + Ln2(D, t)x2(t) +. + Lnn(D, t)xn(t) + Nnn[x(t)] = gn(t) where Lik(D, t) = the linear operators with time-varying coefficients: i=1, 2,...,n k=1, 2,..., n Nik = the nonlinear operators with k= 1, 2,..., n xl(t), x2(t),.., xn(t) = the system output gl(t), g2(t),..., gn(t) = the system input Equation 342 can be put in the matrix form: L11(D, t) L12(D, t).. Lln(D, t) x1 N11 N21'. N x1 gl L21(D, t) L22(D, t)... 2 + 21 2. 2 g Ln(D, t).. Lnn(D, t) x nl Nn x g 78

Institute of Science and Technology The University of Michigan Further, Equations 343 can be put in the following compact form: [L(D, t)]x] + [N]x] = g] (344) where, as defined in Equations 343, [L(D, t)] is a square matrix of the linear operators, [N] is a square matrix containing the nonlinear operators, and x] and g] are two column matrices for the output variables [xi} and the deterministic input functions gi; the independent variable is the time t. If there is only one input g(t) and one output x(t), Equation 344 reduces to the special and simpler case analyzed in Sections 2 and 3. Here, however, we are interested in studying the behavior of the complex system- represented by the matrix Equation 344-when some specified constraints are placed on each of the elements {Lik}, {Fik}, and {gi} for k= 1, 2'...,'n ina systematic manner. The author believes that this will help solve problems of synthesis of complicated systems when deterministic inputs are assumed. To simplify the analysis, the case of two simultaneous differential equation will be considered. (The generalization of the results to systems of n simultaneous equations should follow easily from these arguments.) L11(D, t)xl(t) + N11(x1) + L12(D, t)x2 + N12(x2) = gl(t) (345) L21(D, t)x1(t) + N21(x1) + L22(D, t)x2 + N22(x2) = g2(t) The following eight definitions use the notations presented in Section 2.1: (1) W11(t, u) - the impulse response of the linear term L11(D, t) (2) W22(t, u) - the impulse response of the linear term L22(D, t) (3) F11(t, u, x 1) W11(t, u)N11(X1), and N11(x1) is a function of x1 only (4) F 2(D, t)x2 + N12(x2)], and N12(x2) is a function of x2 only (5) F22(t, u, x2) ~ W22(t, u)N22(x2), and N22(x2) is a function of x2 only (346) (6) F21u, - W22(t, )L21(D, t)x1 + N21(X)], and N2 1(X) is a function of x1 only (6) F21(t, 1' Xl) u=1 \ W22 [211 N2 (7) x10 is the solution of the linear equation L11(D, t)x1o = g1(t) having the same initial conditions as x1(t); that is, x10(t) = W1l(t, u)gl(u) du+terms due to initial conditions (8) x20 is the solution of the linear equation L22(D, t)x1o = g1(t) having the same initial conditions as x2(t); that is, X0(t) = f W(t,u)g2(u)du + terms due to initial conditions By applying the partitioning technique to Equations 345 in such a way that only the linear operators Ll and L22 remain to the left-hand side of these equations, and by then using the above definitions, Equations 345 can be transformed into the following equations in integral form: 79

Institute ot Science and Technology The University of Michigan t t xl(t) = xO(t) - I Fll[t u Xl1(u)]du -0 F12[t' u, x2(u)]du (347) t t x2(t) = x20(t) - j F21[t, U, Xl(U)]du J F22[t, x2(u)]du The relation between circuit configurations and feedback system configurations for the single-degree-of-freedom case has been investigated. Similarly, the n simultaneous nonlinear differential equations can be interpreted as representative of an n-terminal pair (n-port) network, and the equivalent integral equations as representative of a multiloop nonlinear system. To illustrate the latter, consider the system which as two inputs, gl(t) and g2(t), and two outputs, xl(t) and x2(t). The differential equations which govern these functions are given by Equations 345 and can be shown to represent the electric circuit of Figure 24. The circuit shown in Figure 24 has the following equivalences: (1) g1(t) and g2(t) represent input voltages at terminals (1) and (2) (2) x1(t) and x2(t) represent electric current at terminals (1) and (2) (3) [L11(D, t)x1 + N11(x1)] represents the voltage drop on the linear and nonlinear elements of loop (1) due to the flow of current x1(t) (4) [L22(D, t)x2 + N22(x1)] represents the voltage drop on the linear and nonlinear elements of loop (2) due to the flow of current x2(t) (5) [L12(D, t)x2 + N12(x2)] represents, in a general sense, a voltage source in loop (1) that depends only on the current x2 (6) L21(D, t)x1 + N21(x1) represents, in the same way, a voltage source in loop (2) which depends only on the current x! x (t) x (t) Lll and Nl L22 and N 22 - gl(t) L12 and N12 L21 and N21 g2(t) k, I Loop (1) Loop (2) FIGURE 24. A TWO-TERMINAL PAIR, ELECTRIC CIRCUIT REPRESENTATION OF THE DIFFERENTIAL EQUATION 345 80

Institute of Science and Technology The University of Michigan The integral Equations 347 can be shown to represent the multiloop nonlinear feedback systems of Figure 25 or the equivalent conjugate system of Figure 26. The difference between Figures 25 and 26 arises from the difference between the elements in the two minor loops. In Figure 26, the linear time-varying elements of the minor loops, described by W11(t, u) and W22(t, u) in Figure 25, are replaced by the inverse of the nonlinear elements Nll and N22, and are given the names Nll and- N22, respectively. Similarly, the linear operators Lll and Output xl By Input g2(t) ~~~Input g -(t) + N1Output x2 FIGURE 25. A MULTILOOP NONLINEAR FEEDBACK SYSTEM Output xl- 12Input g2(t) -1122 V N12 FIGURE 26. THE CONJUGATE MULTILOOP NONLINEAR FEEDBACK SYSTEM 81

Institute of Science and Technology The University of Michigan L22 in Figure 26 replace the nonlinear elements Nl1 and N22 of Figure 25. When the linear elements are time invariant, the following relations are true: W11(t) = L L(s4 =L 1[Yl1(s)] = 2 — J Y11(s) e+Stds (348) r1) 1T+JOO st W22(t) = L 1 = [Y ()]=2 -1J Y22(s) e ds (349) 122(s ) anJ where L is the inverse Laplace transform operator and Y11(s) and Y22(s) are the linear system transfer functions. Equations similar to 348 and 349 can be derived for time-varying linear elements in terms of the system function [48, 49]. The system function of a time-varying system is defined as the st 31 function Y(s, t) so that Y(s, t) e represents the steady-state response31 of the system to the input g(t) = et. References 48 and 49 show that Y(s, t) and W(t, u) are related to each other by the relations Y(s, t)= W(t, u) e (t )du (350) -00 and U+Jco W(t, u) Y(s, t) es(t-u)ds (351) The salient characteristic of Y(s, t) is that the response x(t) to any input gl(t) may be expressed in terms of y(s, t) and g(t) through the relation x(t) = L- [Y(s, t)G(s)] (352) where G(s) is the Laplace transform of g(t) and, as before, L represents the inverse Laplace transformation. In evaluating the inverse Laplace transform of Y(s, t)G(s), t should be treated as if it were a constant parameter. 31In Reference 49, L. A. Zadeh demonstrates that Y(Jw, t) is the steady-state solution of the differential equation: L(D + Jv; t) Y(Jo, t) = 1. 82

Institute of Science and Technology The University of Michigan 4.3. THE MATHEMATICAL RESTRICTIONS AND THE SYSTEM ANALYSIS (THEOREMS F, G, AND H) In this section, the new admissible class of systems is defined by placing suitable restrictions on the linear, nonlinear, and forcing function terms of the system equations. The physical meanings and the limitations in applications given in Section 2.2 apply here. In particular, the remarks about problems involving the infinite interval of time are assumed to be understood by the reader. To simplify the analysis, the weighted nonlinear functions F11 and F21 are assumed to be functions of x1 alone, and F22 and F12 are assumed to be functions of x2 alone. However, this is not a serious restriction because the systems analyzed here are the same as those given in Figures 25 and 26. The only exception is that the linear branches L21 and L12 are omitted. Now, the following three mathematical restrictions are required: (1) xi0(t) e L2[O, t] O < T < + oo (353) (2) For a given region D; [I x1j ~ dl(t), x2j < d2(t), 0 < u < t S T], the weighted nonlinear functions Fik(t, u, Xk) are assumed to satisfy these conditions: (a) If the two triplets (t, u, Zl) and (t, u, z2) are in D, then ]Fik(t, u, Z1) -k(tk(t, u, z2) ) < (t, t, u ) Z 2 (354) t (b) Fik[t, u, zko(u)] du nik(t) (355) 0 where Kik e L2(0, T), that is, t 2 2 2 0 Kik (t, u) du < aik (t) a (t) and (356) aik2(t) dt < a(t) d(t) < A2< nik < n(t) E L2[0, T], that is ik (t) dt < Cik < n(t) dt C < (357) and i = 1 or 2 (c) It is required that: Xio(t) I < di(t) (358) ]i(t) I< di(t) 83

Institute of Science and Technology The University of Michigan where i = 1 or 2 and where t t X.il(t) = xi0(t) - Fil[t u, x10(u)]du - J F2[t,, x20(u)]du 0 (359) 2o 2(m+2)Am qi(t) Xil(t) l+ c I(t)l 3), m=0 (3) The L2 functions xi0(t), nik(t), and a(t); the time-varying parameters appearing in the system's equations; and the enforced bounds di(t) are continuous and bounded over [0, T]. (360) The method of successive approximations will now be used to prove theorems which will describe some of the properties of the system under consideration. The following inequality will be used in the proof of the theorems: (x + y)2< 2(x2 + y2) (361) Theorem F Let a set of functionsx10, x11,.., Xln} and{x20, x21, X., X2n} be defined by the following recurrence relations: t 4 Xin 1- F11[t, u, xl(u)n-l]du - F12[t, u, x2(U)nl]du (362a) -2n 20 J F2l[t, u, xl(U)n-1]du - F22[t, u, x2(U)nl]dU where n =0, 1, 2,... Then, as specified by the restrictions placed on the system (given by Equations 353 through 360), this set of functions belongs to L2 space and is inside the domain D. Theorem G Under the restrictions specified by Equations 353 through 355, both sequences of iterations {xln} and {x2n } (defined by Equations 356) converge in the mean to the functions yD1(t) and qo2(t), respectively. They also converge in the ordinary sense to the functions &l(t) and z2(t), respectively. Hence, ~l(t)= Pl(t) - xl(t) and (t) 2(t) x2(t) (362b) 84

Institute of Science and Technology The University of Michigan Theorem H Under the restrictions specified by Equations 353 through 355, the functions x1(t) and x2(t) (defined by Equation 357) are the solutions of the system of simultaneous differential Equations 345 and are the unique solutions for a given initial state. Proof of Theorem F Let n = 1 in Equation 356. This gives x11(t) = x10(t) - tJ F11[t, u, x10(u)]du +J F12[t, u, x20(u)]du (362c) x21(t) = x21(t) - F2[t u, x20(u)]du + F22[t, u, x20(u)du 0 The conditions 358 of the restrictions show that the elements x10(t), x1l(t), x20(t), and x21(t) are inside D. However, it must be shown by induction that the two sequences {Xln} and {X2n} are also inside D. Assume that the two elements xInt) and x2n(t) are in D; then x1n+ (t) and x2n+l(t) are defined and are given by t t Xn+l(t) = x10(t) - F 11t, u, Xln(U)]du - { F2[t, u, X2n(u)]du (363):t t X2n+l (t)=x20t) ) = - JF21[t u xin(u) ]du - F21[t, Xln()]du - F22[t u, X2n(u)]du Finally, it must be shown thatxl (t) and x +(t) are in D. Equations 362c give 1n+l 2n+ D E 3lqi i [x11(t) - x1 (t)]2 = tF[t, u, x10(u)]du + F12[t, u, x20(u)]du (364) [x21(t) - x21(t)]2 = F21[t, u, x10(u)]du + F22[t, U, x20(u)]du Making use of the inequalities 361 and the system restrictions 356 and 357 gives [xl(t) - xo(t)]2 < 4n2(t) (365) [x21(t) - x20(t)]2 4n2(t) 85

Institute of Science and Technology The University of Michigan Equations 365 give T [l(t) - x(t)- x0(t)]2dt I Ixll(t) - x10(t) 112 < 4C2 (366) [x21(t) -x21(t)] dt Ix21(t) - X20(t)]2 < 4C2 Subtracting x n(t) from x ln+l(t), as given in Equations 362a and 363, respectively, and then squaring, gives [Xin+l(t)- x = (ln(U)]- Fli[t, u, xlnl(U)]}du J+ {F 12[t' u, X2n(U)] - F12[t, u, X2nil(U)]du < (u lit, xIn(u)] - Fll[t, u, Xlnil(U)] du + 2 F12[t, u, X2n(U)] - F12[t, u, x2nl(u)] du) (367) Making use of Equations 355 and 356 of the restrictions, Equations 367 give [Xin+i(t) 2X(tu)]2 2 du [xn+l(t) - xln(t)] < 2a (t)I [xn(u) -xn- (u)] du+ 2(t) [x2n(u) - x2n-(u)] du (368a) By a similar argument, it can be shown that t t 2n+l(t) -2n(t)] < 2a (t) [2n) - 2n(u)] 2du + 2a (t)f [x2n(u) x2n-l(u)] 2du (368b) Letting n = I in the two equations above and using Equation 366 gives [x12(t) - ll(t)]2 < (22)2 C2a2(t) (369) [x22(t) - x21(t)]2 < (22)2 C2a2(t) 86

Institute of Science and Technology The University of Michigan Again, by following the same argument, it can be shown that t 2 3 22 2 [X13(t) - x11(t)]2 < (22)3 C2a2(t) a2(u) du < (22) C a A (370) [x23(t) - x 2 < (22 3C2a2(t) a2(u) du < (22)3 C a2A 2 Also [xr(t) - x (t)1 2 < (22 )4 C2a (t)A2!x14(t) - x13(t)] - 2! (371) 2 (22)4C2a2(t)A2 [x24(t) - x23(t)] < (2) Ca (t)A2 In general, it follows that xl+2(t) - Xln+l(t) 2 C a (t) (372) 2 22(n+2) A2n C2a2(t) [X2n+2(t) - X2n+l(t)] < n where n = 0, 1, 2,... Now {xt ln+l~ s Jx1(t) + x12(t) -< xl(t) I+ + n+(t) - x(t)n+(t) I and (373) X2n+l(t) I < 1x21(t) + 1x22(t) - x21(t) I +... + X2n+1(t) - X2n(t) 1 Substituting from Equation 372 into Equation 373 gives 00 2(m+2) Am lXln+lxt Ix z I < d1(t) m=0 (374) o 2(m+2) Am m=0 87

Institute of Science and Technology The University of Michigan Thus, xln+l(t) and X2n+l(t) are in D, and it follows by induction that the two sequences {xln(t)} and { x2n(t)} are in D for all n. Since x10(t) and x20(t) are assumed to be L2[0, T], Equations 364 show that xll(t) and x21(t) are L2[0, T]. In general, Equation 372 shows that { Xln(t)} and { x2n(t)} are L2[0, T]; therefore, the proof of the theorem is complete. Proof of Theorem G By the reasoning used to prove Theorem B of Section 2, {xln} and {x2n } can be proved to be two sets of Cauchy sequences or fundamental sequences32. Therefore, from the properties of the L2 space, it can be concluded that two functions oPl(t) and yP2(t) exist in the same L2 space in which {Xln} and {X2n} converge in the mean. In other words, T Limj J [xln (t)(t)] dt= n — oo (375) T2 im Jo [X2n (t) ] dt = 0 n-o 0 Equations 375 show that for every e > 0, there exists a natural number N such that xI xln(t) - eol(t) I I< (376) x2n(t) - sp2(t) < 3E for all n > N. Now, to prove convergence in the ordinary sense, one has to show that the infinite series xll(t) + [x12(t) - xll(t)] + [X13(t) - X12(t)] +... (377) x21(t) + [x22(t) - x21(t)] + [X23(t) - x22(t)] +.. are uniformly convergent. Equations 372 show that Ca(t) 2n+i An-1 [Xln+l(t) - xln(t)] < ()A (378) [x 2n+(t) - x 2n(t)] < C(t) 2n A 32Definition (from Reference 30, p. 169): "A sequence {fn} of points of the space L2 is said to be a Cauchy sequence or a fundamental sequence, if to every E > O there corresponds a natural number N such that I fn - fm lI < c for all m, n > N." 88

Institute of Science and Technology The University of Michigan where a(t), C, and A are finite. One can easily conclude the uniform convergence by using Equations 378 and applying the M test on series 377. From evidence that the n-th partial sum of the infinite series 377 is xln and x2n, it can be concluded that the two sequences { Xln(t)} and {X2n(t)) are convergent. Let Lim x (t)= 1(t) and Lim x2n(t) = 2(t) (379) n-co n n-oo From the use of H. Weyl's lemma (Appendix B), it can be concluded that l(t) = -1 (t) (380) p42(t) = cP2(t) Proof of Theorem H One can easily form the following equations to prove that the limiting functions yi(t) and y2(t) are the solutions of the two simultaneous differential Equations 345: [Xn (t)(t)]2 = {(t) - x10 + tFl1[t u, iln-1(u)]du + F2t, u, x2nl(u)]du + (381) [x2n(t) -y (t)]2 = {2(t) - x20 + F21[t, u, lnl(u)]du + F22[t, u, x2nl(u)]du Then, one can obtain the following equations by integrating both sides of Equation 381 between the limits 0 and T, taking the limits as n - oo, and making use of Equations 375: t t ~l(t) - x10 + Lim S Fl1[t, u, Xln l(u)]du + Lim F12[t, u, X2n-1(U)]du = 0 n-oo n-oo (382) t t qP2(t) - x20 + Lim F21[t, u, lnl(u)]du + Lim 0F22[t,u, 2n(U)] du= 0 n-co n-co Ju Therefore, pol(t) and P2(t) satisfy the integral Equations 347 and are the solutions of the differential equations. 33 Uniqueness of the solution is assured since the sequences {Xln} and {X2nY e L2; however, the proof can be easily demonstrated by considering the system's equations themselves. 33See, for example, Reference 30, page 169, theorem (2). If xln- S~1 and xln- S0{, then Ik|P1 - SD1 |I | S | 1 - XnI + |xlni - qP I || Since the right-hand side has a zero limit, it follows that J]qyl - qo}l = O and =1 = qy{. The same is true for X2n. 89

Institute of Science and Technology The University of Michigan Assume two other solutions po* and yo* and form the following equation: 2 ( O 2 = {It [Fl1(t, u,'P1) - F11(t, u, 1]du t [Fll(t, u, So2) F12(t, u 92)] (du (383) (S~-* = {J [F2l(t, u, 1) - F21((t, u, )]du F[2(t, u, 2,] } Then apply the Lipschitz condition; this gives (-q~ l P1)[2 2 2t 2] < 2a (t) (1 - (p*)2 + ((o2 - (*)j du (384) (-22 < 2 2(t)j [ - ) + (y2 -2)] du Adding Equations 384 together gives 1- 2 +( 2* - 2)] < 4a (t [2 - 9p1) +( 2( 92)2ldu (385) Letting f[(O - 1)2 + (52 - O2)] du = s gives Vp - p1)2 + ( 2 - < 4a2(t)s2 (386) Substituting Equation 386 into Equation 385 and repeating the procedures n times gives JT 22]4n n 2n+2 ( 1O O du < (387) and as n - oo gives T Jo [( ) +2( )2]du 0 (388) therefore, (pl = 9*{ and (p2 = (P*; hence the results. 90

Institute of Science and Technology The University of Michigan 4.4. PHASE SPACE REPRESENTATION AND STABILITY ANALYSIS (THEOREM I) Consider the system of simultaneous, nonlinear differential equations given in Equations 345. Assume that the order of the highest derivative in L11(D, t) is n, and that of L22(D, t) is in. Let X1 = v1 X2 = Vn+ dxl dx2 dt = 2 dt Vn+2 (389) n-1 m-1 dx1 dx2 n-1 n m-1 v dt dt n+m The system of Equations 345 can be transformed to a set of (n + m) first-order equations in the following form: dv dt = f1(vl' v2.. Vn'''' Vn+m' t) dv2 dt = f2(Vl v2' *. Vn'. ** Vn+m' t) (390) dv n+m W dt = f. n+m where {v.}, i = 1, 2,..., n + m, are the state variables which represent the state of the system in the phase space. Each member of the first n components represents the rate of change of the member preceding it. This is also true for the remaining m components. The transformation of the system of Equations 345 to a set of first-order differential equations, as given in Equation 390, canbe illustrated by an example. Consider the following two simultaneous equations: d x dx dx 1 1 2 2 3= 2 + a11(t) —+a- + a12x ++ x+ dt = sint dt (391) ~dx 1x1 dx dx2 dt +a21(t)x +e 1 22+ a22(t) dt + a23(t)x2 = cos t dt2 91

Institute of Science and Technology The "University of Michigan Let n =2 and m =2. Then let x1 = v1 X2 = v3 dxl dx2 dt = V2 dt = V4 Now, Equation 391 can be put into the following form: dv dt = V2 = fl(Vl v2''',' v4) dv2 2 2 dt sin t - v3 - a12(t)v1 - a11(t)V2 = f2(V1. v4) (392a) dv dt = V4 = f3(v1,..., v4) dv4 v1 dt= cos t - - 32 - (= 2f4(v,..., v4) The set of Equations 392a can also be represented in the following matrix form: VI 0 i 0 0 v 0 1 v1 d 2 -a12(t) -a12(t) 0 -1 v2 N2 = + (392b) v3 0 0 1 v3 3 3 v4 -a 21(t) -1 -23(t) -a22(t) 4 -a22( whereN2 = sin t - V3 - v N1= cos t - e In order to investigate stability, as we did for the single-degree-of-freedom case presented in Section 3, we confine ourselves to the free system, represented by Equation 345 in which g1(t) = g2(t) = O. We further assume that the origin of the phase space for the free system is an equilibrium position of uncertain stability. In other words, it is assumed that the function f. (Equation 390) for i= 1, 2,..., n + mis identically zero at v = O andfor all t > t0 = 0. In order to study the properties of the state variables v.(t) of the free system, the system must satisfy a set of conditions which are not necessary, but sufficient and suitable for the present analysis. 92

Institute of Science and Technology The University of Michigan The required system restrictions are (some proofs similar to those presented in Section 3 are omitted to avoid repetition): (1) xi) E L2[0, T] and are assumed bounded over the closed interval [0, T] 0 < t < +oo (393) where (i) = the /i-th derivative of xi (t) with respect to time t, and, 1,., - I for i=, 1,..., m - I for i= 2 (2) In a given domain D [xl S dl(t), Ix2(t) I < d2(t), 0 _< < t < T], the i-th derivative with respect to t of the weighted nonlinear functions, defined as Fik ([t, u, Z(u)], is required to satisfy the following conditions: (a) If the two triplets (t, u, zl) and (t, u, z2) are in D, then Fik( )(t, u, 1) - Fik (t, u, Z2) < Kik (t, u)' z - 32 (394) (b) | Fik( )[t, u, xk0(t)]du nik(t) (395) where Kik (t, u) e L2(0, T), that is, Kik 2(t, u)du < aik 2(t) < a r2(t) L0, T) (396) nk(t) E L2[0, T] The functions aik ar, and nik are continuous and bounded over [0, T] i = I or 2, k = 1 or 2, and /i is as defined above. (3) Condition 358 is still required; that is, I Xio(t) I< di(t) (397) qi(t) Is di(t) where i = 1 or 2 (4) The time-varying parameters and the nonlinear functions |Nik(z) | are assumed bounded over [0, T]. (398) 93

Institute of Science and Technology The University of Michigan The physical meaning of the above restrictions is understood to be the same as in the previous analysis. Now, by differentiating Equations 347 p times and by using the above restrictions and Equations 372, one can proceed, in a manner similar to that used to prove Theorem.E, to find the following equations: l(xn+l (t) - xln )(t)) < 4 2(t) (22(n )CA2(n-1)) (399) (x2n+31 (t) - x() < 4a 2(t) (2 (n 2A)) 2 2 2 where ar is as defined by Equation 395 and A and C are as defined by Equations 356 and 357, respectively. Equation 398 can be used to show that, for the class of systems defined by the above restrictions, the properties of the state variables can be summarized as follows. Theorem I Under the restrictions stated by Equations 393 through 398 for the system of two simultaneous, nonlinear differential equations given in Equations 345, the following properties of the state variables v.(t), i = 1, 2,..., n + m, are true: (1) vi(t) E L2[0, T] (2) vi(t) are bounded by corresponding L2[0, T] functions ji(t), as given by x (i-(t)I+ In1 (t)I + nl2 (t) + mlar(t) for i= 1, 2,..., n pt(t)= (400) x20 (i21 (t) t+ n+22 n(t)I + In+mlr(t)I for i=n+l,., n +m 1-n-1 i-n4 where m = 4C ( 2A1 + (a)n (401) and C, A, and ar(t) are as defined before. (3) vi(t) are uniformly continuous over [0, T] (4) If T - +co, then it follows that Limv.(t) = 0 (402) t-co i 94

Institute of Science and Technology The University of Michigan From Theorem I it can be concluded that at any time t for which 0 <t < T, the following inequality is true: vi(t) < pi(t) i = 1, 2,..., n + m (403) where /i(t) are as defined in Equations 400. But v.(t) are the limit of sequences {vin(t)} where i = 1,..., n + m and n = 1, 2,..., oo; which belong to L2[0, T]; therefore, we deduce that the L2 norms Ilvi(t) I are bounded. This can be concluded from the following corollary: Corollary: The norms of the elements of a convergent sequence in L2 are bounded [30]. However, for the class of system under consideration, the upper-bound functions i(t) can readily give the required upper-bound norms of the state variables so that: Ivi(t) I I< I13i(t) I i = 1, 2,..., n + m (404) However, it should be noted that where P(t) = {ji(t)}, the upper-bound norms I [3(t) I I are functions of the initial state V0 = {vi0} at t = 0 and i = 1, 2,..., n + m, so that an exact definition of boundedness can be stated as follows. Definition [31, p. 388]: A motion is said to be bounded for every (V0, t0) if there is some constant P(_0, to) so that I Jv(t; v0, t0) < ~3 for all t > t0. The motion is equibounded if 3(v0, to) < 3(r, t0) for all ~I v0 I I < r. Now, by using the same two definitions of the norms that were used in Section 3: n+m 2 v. it can be concluded, just as it was for the single degree-of-freedom case, that the class of nonlinear systems considered here is asymptotically stable. 4.5. EXAMPLE CALCULATIONS As an illustration of the method presented in this section, consider the multiloop system shown in Figure 27. In this system there is only one input gl1(t) and two outputs xl(t) and x2(t). The two minor loops are coupled together by two nonlinear elements, each of which has the 12 form -x. The boxes marked "b" (b > 0) represent variable gain amplifiers with amplifier 95

Institute of Science and Technology The IUniversity of Michigan x (t) FIGURE 27. A MULTILOOP NONLINEAR FEEDBACK SYSTEM WITH FOUR NONinearLINEAR ELEMENTS, x2 t gains b. Th l e of the amplifier gains b, so that the outputs x1(t) and and 28b Noalso equivalent tonlinear Nonlinear two-terminal pair network of Figure 29. onlinear differential equatiNonlinearons: \ dt + +/ + 2) =be2t (406) + X + -2) 0 dx1 -2t b 2 b 2 dt +x1 = be - x1 -x2 (407) 2b 2 b 2 dt +X2 = 2+ 2X1 - X2 96

Institute of Science and Technology The University of Michigan xl(t) -2t x2(t). =e R = 1 L = 1 R =i R L =1 R =1 b N - Nonlinear- Nonlinear 1/V~ S n21 2 L I v.=2C1 2 2 2 FIGURE 28. A TWO-ULTILOOP NONLINEARCIRCUIT EQUIVALENT TO THAT OF FIGURE 27 L 2l n 2' n 12 2L= ii2 v2(t)= be 0 I (t) COUPLING 97

Institute of Science and Technology The Oniversity of Michigan The eight definitions introduced in Section 4.2 are evaluated for this particular system as follows: Wll(t, U) = e(tu) W22(t, u) = e(tu) x10(t) = e + b(e -t 2t [at= 0 x10(0+) = 1] x20(t) = e, [at t = Ox20(0+) = 1] (408) Fll[t, u, X1(U)] = 2e (tU)x1 (u) F12[t, u, x2(u)] = b2 e(tu)x22(u) b -(t-u) 2 F22[t, U, x2(u)] = e x2 (u) b -(t-u) 2u) F21[t, u, Xl(U)] = e x 1 (u) If we choose a domain D [IxIl < 1, |x2| I 1, O <u<t Sl], condition 354 will give K.. = b e(tu) for i = or 2 for j = 1 or 2 (409) aij = Kij2(t, u)du= (1- e a2t) = 2(t) (410) t 2 -2e A. 2 = 2(t) dtb2 e =0.283b2 (411) Thus A..= A = 0.53b IJ Applying conditions 355 and 357 gives n21(t) = nl1(t) = F11[t, u, x10(u)]du (1 +b+ ) e-t _ (1 + b)2 e2t +b(1 +b) e3t -4 (412) 98

Institute of Science and Technology The University of Michigan b (,- e (413) 12(t)= J F22[t u, x20(u)]d = (et - e2t (413) C2 = Cll n (t)dt - (1 + b) 0.432 - (2/3)( (1 + b) + (414) C22 = C12= n22 (t) dt = 0.004b2 (415) Since for all t 0 < t < 1, xl0(t) I should be less than or equal to one, b can be chosen less than one. Thus, b can be chosen to satisfy the inequality 0 < b < 1. Take b = 0.5 as a representative 2 2 2 value. For this particular value of b, it is possible to choose C = C22 = 0.004b = 0.0002 or C = 0.014. It also follows that 2A = 0.53. A check of inequalities 358 shows that they are satisfied. Now, use Equation 400 and notice that an upper-bound function for x1(t) is obtained for i = 1, and an upper-bound function for x2(t) is obtained for i = n + 1 = 2. It can be shown that for the system under consideration and in the given interval, the required solutions x1(t) and x2(t) can never exceed the upper-bound functions Pl(t) and:2(t) as given by Ixl(t) I-< l(t) = x10(t) I + nll(t) I + In12(t) + mIal(t) (416) x2(t) I 32(t) = X20(t) I + In21(t) I + in22(t) I+ m|a2(t) I where m is given by m =4C Z(2A)k (417) k=0 In our case, m = 0.003 and xl(t) < 0l(t) = 2e-t - e -2t + 0.01(1 - e2te1/2 (418) x2(t) | 2(t)= e-t -t -2t) -2t 1/2 1x2(t)12(t) = e +.5(e - e +0.1(1-e ) for 0 <t < 1. 99

Institute of Science and Technology The University of Michigan In order to check the validity of inequalitiies 418, it is necessary to find the exact solutions x1(t) and x2(t) for the system under consideration.34 The method of analysis used in this section readily leads to an approximation of the exact solutions of the system to any degree of accuracy by use of numerical techniques programmed for a digital computer. Thus, by substituting from definitions 408 into Equations 362a, the following recurrent relations can be obtained: x1nl(t) = e t + 0.5(e - e 2 - et e xin (U) + X2n (u du (419) -t t X n+l(t)=e + ex (u - X22(u)]du Since Theorem H proves that the recurrent relations 419 will yield the exact and unique solution as n - oo, the convergence of the process is guaranteed. However, by stopping at the n-th iteration, it is easy to prove that the error 6, in the mean square sense, for the given class of multiloop system can, in general, be given by the following equation: n-= C2A)= L (/)k] (420) where n2= (X - xn)2dt (421) 2 2 and C and A are as defined before. The derivation of Equation 420 follows essentially the same procedure used to obtain Equation 170. Substituting for C2 and A2 in Equation 420 gives 62 0.00063(053) (422) n 0. (n - 1)! and for n = 1 6 2 = 0.00017 (423) 34 The inequalities can easily be proven true because, by inspection, Equations 406 for b = 0.005 has solutions xl(t) = x2(t) = et. 100

Institute of Science and Technology The University of Michigan Thus the solutions can be approximately given by Xll(t) and x22(t) as follows: xl(t) -Xll(t) = x10 - n11(t) - n12(t) (424) x2(t) x21(t)= x20 + n21(t) - n22(t) Substituting for the quantities in Equation 424 gives -t -2t -3t -4t x1(t) x11(t) = 0.85e + 0.318e - 0.189e + 0.021e (425) x2(t) -x21(t) = 1.15et - 0.318e + 0.189e - 0.021e-4t for 0 t < 1. Note that the approximate solution 425 is valid only for the finite interval 0 < t < 1. Thus, the problem cannot be worked for the infinite interval of time (for the same reasons that apply to example one, Section 2). If an analysis of the problem is required for the infinite interval of time, some modifications are necessary. For example, d (t) and d (t) can be chosen so that -at -at 1 2 d1(t) = e and d2(t) = e 2, where a1 and a2 are finite positive numbers. This has been shown in detail in example three, Section 2. 4.6. CONCLUSIONS AND REMARKS The analysis developed in Sections 2 and 3 was extended to include systems whose dynamic behavior can be represented by a set of simultaneous nonlinear differential equations that are difficult to reduce to the canonic form of a single-degree-of-freedom nonlinear equation. A class of multiloop systems whose exact analytic solutions can be obtained to any degree of accuracy within a prescribed error was found to exist. The type of error used to estimate the system response was the mean square error, since it was best suited to the method of analysis. The mean square error can be used when the instantaneous value of the error in the approximate system response is not as important as, for example, the value of the integrated square error in the interval of interest. By placing suitable restrictions on the different terms appearing in the system equations, it was found that the state variables [vi(t), i = 1, 2,..., n + m] of the system belong to an L2 space, and that each of these state variables is the limit of a sequence of iterates which also belongs to L2 space. Convergence of these iterates within the given restrictions was satisfied in two cases: (1) the mean and (2) the ordinary sense. Note, however, that convergence in the mean of a sequence of functions which belong to L2 space does not imply convergence in the ordinary sense. The converse is also true. (This fact is made clear by examples in Appendix I.) 101

Institute of Science and Technology The University of Michigan Because the norms of a convergent series in L2 are bounded, one can be certain that a system which operates within the given conditions will have a response which is bounded. If the interval of interest T is taken as infinity, and if the system satisfies the given conditions in that interval, the system will return to the equilibrium position-assumed to be the origin of the phase space-and give rise to asymptotically stable systems, in the sense of Lyapunov. As in the case of the single-degree-of- freedom system analyzed in Section 3, two definitions of the norm of the response functions were used to investigate the behavior of the system response at any instant and to estimate, on the average, the upper-bound functions that belong to L2 space. As a by-product of this work, some interesting relations between n-terminal pair circuit configurations and multiloop feedback system configurations were developed and fully explained. Also, the output responses of the feedback systems were shown to be equivalent to the responses of the currents in the circuits. Finally, an example was given to illustrate the method of the analysis presented. SYSTEM ANALYSIS BY NUMERICAL TECHNIQUES The preceding three sections presented an analysis of a class of physical nonlinear systems represented by either a single-degree-of-freedom equation or a set of simultaneous, nonlinear differential equations. Attention was focused upon the properties of the system response when certain conditions were satisfied and upon the problem of how to choose system parameters to yield certain requirements. These considerations are an important step toward solving problems of synthesis. The exact analytical solution for the class of systems under study was shown to be obtainable from Picard iterates, and the convergence of these iterates was assured. After a finite number of iterations, say N, is obtained from an appropriate formula, then the N-th iterates will provide the solution for a predetermined degree of accuracy. This section discusses the problem of obtaining a numerical solution for the class of systems under consideration, by using a digital computer to numerically evaluate the iterations. This reduces the complexity of a nonlinear system calculation almost to that of numerical integration. However, the accuracy of the method depends on the accuracy with which the, superposition integrals appearing in the iterates are evaluated. Therefore, the problem now is to obtain approximations for integrals, derivatives, and intermediate values (interpolations) of the continuous function x(t), at any time t, when given only a table of values xs at discrete stationary points t = ts. This necessitates a discussion of interpolation, numerical integration, and differentiation. However, discussion here treats only the situation in which the spacing h between station points (h = ts+l - ts) is uniform, 102

Institute of Science and Technology The University of Michigan i.e., h is constant and independent of s. Although there are various methods for solving ordinary differential equations numerically, a comprehensive treatment of each is beyond the scope of this work. Therefore, only the most pertinent are discussed; in particular, the fourth-order Runge-Kutta method will be used. This section also discusses theoretical considerations of the accuracy of the iterative procedure and a method that has been devised for estimating "inherited" errors in terms of "truncation errors." This, of course, enables one to choose an interval h so that a predetermined accuracy can be maintained. 5.1. INTRODUCTION By following Albert A. Bennett's classification [6, p. 65], approximate methods of integration of differential equations can be separated into two main types. Type 1 The desired region of definition is assigned in advance or made as large as possible. For this region an assumed type of solution with infinitely many parameters is progressively evolved, until a satisfactory preassigned degree of agreement is assured. This method includes, for example, power-series solutions by undetermined coefficients, if convergence through the region can be proved or is assumed. It also includes the Picard method of successive substitutions. Type 2 An approximate solution, which is usually a polynomial of not more than sixth order, is chosen in a form suitable for each interval of not too great extent. The parameters in the solution and the length of an interval starting at the initial point are chosen so that a satisfactory preassigned degree of agreement is assumed in this interval; then, this procedure is repeated for each interval until the desired preassigned total region in the real axis is covered. Such a method can be called a step-by-step method, in the restricted sense. The use of quadrature formulas, expressed in terms of differences, to proceed from one step at ts to the next at ts+l is essentially a type 2 method. Adam's method, Milne's method, and Multon's method are all typical examples [27]. One other interesting method that belongs to type 2 but differs from the processes which use quadrature formulas is the Runge-Kutta method [13]. This section deals principally with the type 1 method: the use of iterative procedures for the systematic approximation of the partition method used to analyze the class of nonlinear systems under study. A comparison with similar methods of analysis developed by Pipes [32], Stout [36], and Ku, Wolf, and Dietz [42] for the approximate solution of similar nonlinear integral equations of the convolution type is given. The accuracy of most of the above methods depends chiefly on the accuracy of the numerical evaluation of certain integrals and derivatives. (A brief discussion of some of the problems involved in handling discrete functions is presented 103

Institute of Science and Technology The University of Michigan in Appendix J and, a more detailed discussion, in References 9, 27, and 34.) Some of the errors that can occur in the analysis are truncation errors, round-off errors, and discretization errors. Round-off and discretization errors are defined and discussed in Section 5.4 and truncation errors are defined in Appendix J. 5.2. THE ESTIMATION OF TRUNCATION ERRORS BY TWO MEASUREMENTS It is not always possible for a computer to have an explicit formula which gives the value of the truncation error committed in the evaluation of a given result. Since truncation errors are important in providing tentative preliminary estimates of actual error propagation, a fairly simple method is given for the determination of truncation errors and the determination of the exact value of a given quantity in the computation. The method is applicable to the computation of definite integrals and to other computations as well. The Method Let a quantity I be approximated by a quantity R1 with a truncation error E1, and by a quantity R2 with a truncation error E2, in such a way that E2 =K E1 (426) where K is any given constant number. *Then, R1 - KR 1 - K (427) and R1 - R2 1 1 (428) E1 = 1 1 K The Proof The proof of Equations 427 and 428 proceeds as follows. I=R1 +E1 =R2 +E2 (429) Using Equations 426 and 429 gives (I - R1) = K(I - R2) or I(1 - K) = R1 - KR2 104

Institute of Science and Technology The University of Michigan which proves Equation 427. Also, E1 =I-R1 R1 - KR2 1R 1-K -R1 R1 - KR2 - R1+ KR1 R1 R2 1 1 -K 1 1K which proves Equation 428. Equations 427 and 428 will now be used to obtain information about the truncation errors in the following three cases. (1) Application to the Trapezoidal Rule Let a definite integral I be estimated by the two quantities R1 and R2 by using the trapezoidal rule and two different values of h in each computation. Equation 547 of Appendix J shows that the truncation error is of the order of h; hence, h2 = 1/2 hi, corresponds to k = 1/4 and, by its substitution into Equations 427 and 428, gives R1 - R2/4 4R1 - R i= 1 2/4 1 2(430) R -R R - R2 E 1 2 1 2 (431) 1 4-1- 3 (2) Application to Simpson's Rule Equation 547 of Appendix J shows that the truncation error in Simpson's rule is of the order of h. Therefore, two measurements with h2 = 1/2 h1 give a value of k = 1/16. Now, substitution of this value into Equation 427 and 428 gives R - R 16R -R 1 16 2 1 2 (432) I=,, I (432) 1' 15 16 R1 - R2 R1 - R2 1 2 16-1 15 Equation 433 is in agreement with the result obtained in Reference 35 (p. 175). 105

Institute of Science and Technology The University of Michigan (3) Application to the Runge-Kutta Process The Runge-Kutta method of the fourth order' uses the first five terms of the Taylor series and, hence, includes an error of the order of h. If two measurements of the ordinate x1 (as defined in Figure 30 are obtained, first in one step through an interval h2 and second in two steps through two intervals, each of length hl, so that h2 = 2hl, the result will be a value of k = 1/16. Therefore, substitution of this value into Equations 427 and 428 gives: 1 RI -- R2 16R1 - R I = 1 1 (434) 1 1 15 R1 -R2 R1 - R2 E1 K 1 1 (435) where R1 and R2 are the two computed values of X1. Equation 435 agrees with the result obtained by Ince in Reference 13 (Appendix B). This equation was used to solve a nonlinear problem on the digital computer by means of the Runge-Kutta process, and its use was helpful in obtaining satisfactory results as well as in providing a means for checking the computed results during the computation process. Since the formulas, methods, and definitions usually associated with numerical data have now been explained, the rest of this section will present the application of these rules to the numerical solution of problems pertinent to this report. x1 x1 t t (a) (b) FIGURE 30. THE MEASUREMENT OF x1. (a) Two steps each of interval length hl. (b) A single step of interval length h2 = 2hl. 5.3. THE SYSTEMATIC APPROXIMATION TO THE PARTITION METHOD First, consider the nonlinear system defined by Equation 96 together with the restrictions 102 through 109. For this class of systems, the solution was shown to be given by the following 2 3 4 =n+ X +hx' x ++ x'" + x +O(h) n+ n 2 n 3! n 4 n 106

Institute of Science and Technology The University of Michigan iterative procedure: Xm+1 (t) = x0(t) - W(t, u) N(xm (), xm(),, )u (436) m'Xm,xm du (436) where m = 0, 1,2,.., p,...; O < t < T For a given mean square error q, the p-th iterate was shown to approximate the system solution in the mean square sense where the value of p, when N is a function of x alone, is obtained by solving for m as follows: 2 C 2m A (437) k=O 2 2 and C and A have the meanings explained in Section 2. To calculate the integrand of Equation 436 for each step of iteration, further expressions are needed for' Xm'..,andx (n-1) are needed for x, xm,.., and xm, where n is the order of the system equation. Therefore, Equation 436 can be replaced by its equivalent relation xm+l (t) = x (t)- W(t U) N(xm xm, u) du (438) i = 0, 1,..., < n- 1 as was shown in Section 3 for a similar equation. The solution of Equation 96 requires the evaluation of x(i)(t) and the substitution of these functions in Equation 438 to obtain a first approximation for the solution x1(i) and so on for further approximations. This process is greatly expedited by using a digital computer. Therefore, the interval in which the solution x(t) is desired is divided into continuous disjoint subintervals of equal length h. The solution, then, has the following approximate values: x (t), xm(1)(t), x0(t), x0(1)(t), etc., for t = tk = kh; these are denoted by xmk, xmk( Xk, X0k( ) etc. In particular, W(t, u) and N(x, xm, x(2),..., t) are defined as follows: W(t, u) = Wks, where (O < s < k < (439) N(t) = N [xm(t), xm()(t), ) xm( )(t), * = t* * t = Nm (440) 107

Institute of Science and Technology The University of Michigan Then, for t = tk, Equation 438 becomes m+(i) (tk) = (i)(tk) - i w(i)(tk, u) N (u)du (441) If the superposition integral appearing in Equation 441 is approximated by any particular rule of numerical integration (the general quadrature formulas are explained in Appendix J), Equation 441 can be written approximately as follows: x )(i) xOk(i) - h y W (i)N (442) s=O where ys is a constant determined by the particular integration rule employed. (See Equations 541 through 543 for the special cases of the trapezoidal rule, Simpson's rule, and Weddle's rule.) For a more convenient representation, Equation 442 is give n in the matrix form: X(m+l)k x Ok Wko Wkl.' Wkk Y0 Nm0 (m+l)k = Ok - h kO Wk kk 1 Nml (443a) (n-l) xn-1) (n-1) x(m+l)k Ok k(n 1) Wkk W k Nmk where each term is as defined before. Equation 443a gives the (m + 1)-th iterate in terms of the previous and known iterate of order m and gives the functional values for each discrete value of K, (k = 0, 1,..., T/h), where T is the period of interest. For example, the first iteration can be obtained in terms of the known discrete values xk xOk ).., xk(n-l) in the following first-order iterative matrix: xlk Ok Wk Wk... kk YO N00 x(1) x (1) h W k(1) W (1).W (1) Y1 N (443b) xk(n- (n-i) (n-i) (n-i) ( Yk N 108Ok kO ki kk 71 1 108

Institute of Science and Technology The University of Michigan Then, from this matrix, a second-order iterative matrix can be obtained: X2k XOk Wk0 Wkl... Wkk N (1) X( 1)Ok W W Wkk - h k k kk 1 (443c) (n- (n-i) 1 (n-i) W (n-i) 1 (n-i) Yk Nk) X2k- _O0k n kO k1 W kk and so on, until the p-th order matrix which closely approximates the required solution to a given degree of accuracy is obtained. Obviously, the column matrix [x0k(i) and the square matrix [Wks] appearing in the right-hand side of Equation 436 remain the same in each iteration. 5.4. THE ACCURACY OF THE APPROXIMATE SOLUTION The accuracy of the method described above depends largely on the particular method of integration adopted and on the particular problem being solved. The errors in the quantities mk, which are the i-th derivatiative of x evaluated at the time tk in the m-th iteration, are defined by = x (i) - x () (kh) (444) (i,k) mk m These errors are due to the truncation errors in the approximate evaluation of different integrals in the k-th step and the "inherited" errors from similar approximations in the previous steps. The error defined in Equation 444 is sometimes called the discretization error [9, p. 16], a name that was used by Crandall. The discretization error describes the actual propagation of the error at each discrete time in the calculation. Frequently, the discretization error is difficult to evaluate; consequently, the truncation error is sometimes used as a crude estimate for the discretization error. A qualitative argument [9, p. 170] is used for its justification (for example, see Reference 36). Generally, it can be said that at the end of n steps of calculations the total error should be of the order of n times the truncation error of a single step [9, p. 167]. It is known that, for a given fixed total interval, the number of steps n is inversely proportional to h. This will reduce the order of dependence of the total error at the end of the interval on 36 T the total interval of interest ActuallY, n h the subinterval width 109

Institute of Science and Technology The University of Michigan the size of the increment h. The following example illustrates this. Assume the use of the simplest form of Euler methods to obtain the numerical integration of the first-order differendx tial equation dt = f(x, t). For a given reasonable value of h, we have tial equation 2 xn+1 x nh(h + O(h) (445) n+ln \dt=nh where O(h ) means that the added terms are of the order of h. From Equation 445 it can be seen that the truncation error in each step of the calculation is of the order of h2. At the end of T h subintervals, the error will be h 2 eKt. = h..() (446) which is of the order of h only. An exact estimate of the discretization error en [9, p. 168] was found to be n p (+ a epT -1 - 2 (pTePT) + O(h2p2)] (447) n 2p at 2 avg. where p =-. ax' Examination of this equation shows that if the value of hp is small, the discretization error n is proportional to h, and this agrees with the previous intuitive argument. Now, we can conclude that the discretization error can be decreased by using smaller increments of h, but the number of steps and hence the amount of computation required to cover a given interval are increased. When a large number of steps is required, there is a danger that round-off errors37 will build up to substantial proportions. However, since the magnitude of the round-off errors can be controlled by the number of decimal places to which the computation is carried, we can consider round-off errors negligible in comparison with the other errors mentioned above. Note that in the last example, and in the special case in which the maximum possible round-off error r at the end of a fixed interval is inversely proportional to h and the maximum of the discretization error n is directly proportional to h (Figure 31), the minimum of the maximum of the total error E = + en occurs when the round-off error E and the truncation error ET of a single step are of equal magnitude; that is, when [9, p. 169] E = E Note that E equals r T' r the greatest possible round-off error at each step and ET equals the greatest truncation error at each step. Thus when the round-off errors are neglected, the truncation errors are important in the evaluation of the discretization errors ~(ik)m' as defined by Equation 444. There37Round-off errors are errors committed in rounding off results to a finite number of digits. 110

Institute of Science and Technology The University of Michigan Total Error Curve &(at t = h) Max e T h Optimum Step Size h =n FIGURE 31. A PLOT OF 8r AND En VERSUS h fore, the maximum absolute value of the truncation error for each step in the solution must be much less than the actual maximum discretization error in the solution. However, a preliminary value of h can be derived from an estimate of the truncation error. In the demonstration that follows, the trapezoidal rule of integration is used. Although it is less accurate than other methods, it has been found to yield satisfactory results [36, 42]. Equation 545 is used for the truncation error evaluation and in each iteration there are integrands z of the form: z = w(i)(t, u) Nm(u), where i = d (448) If one takes the first derivative of Equation 448 with respect to the variable u, dNm (u) d(i) d u ( W i)(tu) du (449) lets t = nh, substitutes u = nh and u = 0 in Equation 449, and uses Equation 545, one obtains the following expression for the truncation error: ET[xm+(i) (nh h ()(nh, nh) N(nh) + Nm(nh) [W(i)(nh, mh)]' (450) - W(i)(nh, 0) Nm(0) - N'm(O) [W(i)(nh, O)]} where the primes mean differentiation with respect to u. (Note that the truncation error resultb from an approximation of the superposition integrals appearing in Equation 441 when the (m + 1)-th iterate is evaluated using the trapezoidal rule.) For the special condition in which we have a nonlinearity N which is a function of x only (therefore, only the iterates of x, and not its deriva111

Institute of Science and Technology The Ulniversity of Michigan tives, are required) and a linear part with a constant coefficient (therefore, W[t, u] = W[t - u), Equation 450 reduces to: 2 E[xm+1(nh)] <2 [W(O)Nm(nh) - N' (nh)W'(O) - W(nh)N' (0) + Nm(0)W'(nh)] (451) Equation 451 is in agreement with the one derived by T. M. Stout in Reference 36. Of course, the calculation of the derivatives of Nm(u) necessary for the evaluation of Equations 450 and 451 can be obtained to any degree of accuracy by the methods described in Appendix J. However, a satisfactory estimate of the derivatives can be obtained from the following equation: N (n + l)h - N (n - l)h 2 N' (nh) = m 2h (452) M 2h +0(h) (452) 2 if the necessary values are available. The error in Equation 452 is of the order of h Now that the truncation errors in the approximations are known, the propagation of the discretization errors as given previously by Equation 444 can be studied. The truncation errors ET[Xm (i)(nh)] are denoted by E(ik)m, so that the exact (m + 1)-th iterates of the differential Equation 96 at t = kh can be obtained by applying the following equation: Xm+(i) (kh) = x0(i)(kh) - h Wks(i)N (sh)+ E(ik) (453) s=0 Substituting Equations 453 and 442 into Equation 444 gives (ik) ml 0k() x( ) (kh)1 - h4 Y/sWks[Nms - Nm(h) - E (i,k)m where, as before, the y's depend on the particular rule of integration used and i < n = the order of system. Equation 454 gives a recurrent relation that represents the propagation of the errors in the discrete values of any of the iterates, and that greatly controls the accuracy of the method presented. The method is applicable when the forcing.function g(t) is specified graphically, tabularly, or by an explicit mathematical function. However, when g(t) is given by an explicit mathematical function, the accuracy of the method can be improved by using the exact response of the linear partitioned terms x0( )(kh), so that [x (i) -x0(i)(kh)l= (455) 112

Institute of Science and Technology The University of Michigan and Equation 454 reduces to = -hi Yswks[Nms (sh)] - E(i k)m(456) The errors in the discrete values of the nonlinear function represented by the terms [Nms - Nm(sh)] are nonlinear functions of the errors C(i,k)m' Equation 456 can be simplified if the approximate discrete values of the nonlinear functions N are expanded in a Taylor series m s about the exact values of x(t) and its derivatives appearing in the function. Therefore, we have Nmk= N (kh) + Aj (Jk)(457) in which A = [W;] (458) Jm ax( (t) t=kh so that [Nmk- Nm(kh)] = AJm(J,k) (459) 04) Nm(t). The quantities AJm usually are functions of x (J)(kh); but constant values for AJ can be * m m J assigned over a limited range of the expansion. When this is done, Equation 456 reduces to m+i S=O -i MM (i k) = -h i LE s; AJ 6 ~k(J.s)j iS - E (460) where m + 1 = 1, 2, 3,.., q, (q is the q-th iterate) i = 0, 1,..., n - 1, (n is the order of the system) k = 0, 1,... T (T is the total interval) Again, (ik)m+l denotes the discretization error in the (m + 1)-th iterate of the i-th derivative of the solution x(t), at the k-th discrete station point. The truncation error E(i,k)m varies, of course, at each point tk = kh for each value of m, but a further simplification can be obtained by assigning a value Ei(kh) for each of the quantities E(ik) for all m. In this case Equation 460 takes the more simple form ~(i)m+1 h Y= r A~ ~(Js)m] W(ks)(} -Eik (461) 113

Institute of Science and Technology The U niversity of Michigan Changing the order of summation in Equation 461 and approximating the sum over s by an intergration process reduces this equation to the new form E. (t) = [-Ei(t)] - Aj j W( (t, u)J (u)duj (462) 1(m+l) Letting AJW(i)(t, u) = Wi(t, u) (463) reduces Equation 463 to the form J<n-1 t i(t)(m+1) =[ -Ei(t)] - WiJ(t, u) UJ (u)du (464) J=0 It was proved in Section 4 that at the limit at which n goes to infinity, Equation 464 yields the following equation: J<n-1 t (t) [-E(t)] W. (t, u) FJ(u)du (465) 1 J0 wherei = 0 1, 2,..., < n- 1, < t< T Equation 465 is know as a linear system of Volterra's equation of the second type [24, p. 17]. Its solution represents the error propagation at the end of a sufficiently large number of iterations. For simplicity, the solutions can be assumed to be adequate representations of the discretization error propagation at the end of the required number of p iterations necessary to approximate the system solution by the method presented. For the special condition in which the nonlinearity N is a function of x only, Equation 465 reduces to ~(t) - X W(t, u) &(u)du = -E(t), 0 < t < T (466) 0 where X -A0 Aax(t) (467) Equation 466 is a linear Volterra equation of the second type,38 the unique solution of which can be given at once by the formula: 8(t) - [-E(t)] + H(t, u, X) E(u)du (468) 38If E(t) = 0, Equation 466 reduces to an equation of the first type [38, p. 10]. For the purposes of this report, it is assumed that Ei(t) will always have a certain value. 114

Institute of Science and Technology The University of Michigan where the "resolvent kernel" H(t, u, X) is given by the series of iterated kernels as follows: oO -H(t, u,:) = X Wn+(t, u) (469) n=0 Series 469 converges almost everywhere. The resolvent kernel is known to satisfy the integral equation: W(t, u) + H(t, u, x) = X X W(t, z) H(z, u, X)dz = X 5 H(t, z, X) W(z, u) dz (470) Note that the iterated kernels W (t, u) are defined as t W +( (t, u) =j W n(z, u)dz n1, 2, 3... (471) where W(t, u) = W(t, u) (472) From the error Equation 468, it can be concluded that the error in the solution by the procedure presented depends primarily on the truncation error E(t), and secondarily on the nature of the particular problem in the region of interest. The dependence on the region of interest is represented by the X term39 in Equation 468. If the truncation errors are reduced significantly -and this can be done by using more accurate rules of integration-then the effect of the X term will be negligible. For the general case, Equation 461, 464, or 465 can be used to give information about the error propagation. The choice of the particular equation depends on what is required: accuracy or simplicity of analysis. 5.5. AN EXAMPLE CALCULATION USING THE DIGITAL COMPUTER The use of the numerical analysis technique given in Section 5 can be illustrated by the following equation, in which the technique is applied to a nonlinear system (given in Section 2): dx 2 -2t dt + x + x = e where x = -1 at t = 0 (473) 39 See Equation 467. 115

Institute of Science and Technology The LUniversity of. Michigan Partitioning this equation at the x term and using previous notations and definitions, adopted in Section 2, results in.the following recurrent relation: Xn+l(t) =x0(t) -X e(tu) x2 (u)du (474) where -2t x0(t) = -e (475) Now, the values given in Section 2 for C (0.025) and A (1.14) can be substituted into Equation 170 to produce the following equation: 2 (0.025)(9.6)(1.14)n q (n - 1)] (476) Equation 476 is used to construct Table IV, which gives the required number of iterations n for 2 a given mean square error q. Table IV shows that it is reasonable to perform the first seven iterations of Equation 474 on the digital computer. If one assumes that the trapezoidal rule is used to evaluate the different integrals, Equation 451 can be used to give an estimate of the truncation error committed in evaluating the integral appearing in the first iterate at the end of the interval of interest T, which is assumed to be unity. For this purpose the following quantities are needed: W(t) -= et N(nh) =x02(t) = e4t (477) -t N'(h) -4e4t W'(t) = -e N(nh) = e TABLE IV. THE MEAN SQUARE ERROR q2 FOR DIFFERENT VALUES OF ITERATIONS n n = the number of q2 = the mean iterations square error 2 0.31 3 0.118 4 0.0635 5 0.018 6 0.004 7 0.00075 116

Institute of Science and Technology The University of Michigan Evaluation of these quantities at t = 0 and t = 1 gives W(O) = 1 W(1) = 1/e W'(0) = -1 W'(1) = -1/e (478) N1(0) = 1 N1(1) = l/e N'1(0) = -4 N' 1(0) = -4/e4 Choosing h = 0.05 and substituting this together with the quantities of Equation 478 into Equation 451 gives El[x(l)] < 0.000218 ('479) After the required number of iterations, an adequate estimate of the propagation of errors can be obtained by considering Equation 466 and assigning reasonable values for E(t) and X. Thus, 2 N(x) = x (480) and -A0 x= A ax -2x(t) (481) Although X varies with time, a constant value is assigned over the limited range of the solution. Also, E(t) = E(kh) varies at each point in the solution so that E(0) = 0 (since the initial point is known exactly from the initial conditions). Thus, a function E(t) can be obtained if the truncation error is known at each discrete point, in the first iteration, for example, and is considered as the required function E(t) over the required number of iterates to be substituted in the error Equation 466. However, the upper bound of the error &(t) can be obtained quickly by assigning to E(t) a constant value equal to its worst value, such as that given by Equation 479. Therefore, the actual error in the solution of the given example will be bounded by the solution of the following equation (obtained from Equation 466), after the appropriate substitutions for W(t, u) and E (t): ~(t) + A 0 e-(t) E(u) = -0.000218 ('482) where X is a constant parameter. Equation 482 can easily be transformed to a linear differential equation. Let Et (t) = 6(t) et e (483) f1(t) = -0.000218 = -0.000218 e+t -t1 117

Institute of Science and Technology The University of Michigan Dividing both sides of Equation 482 by e-t and making use of Equations 483 gives l(t) + A0 (u) du -0.000218 e(t (484) Differentiating both sides of Equation 484 gives dl(t) +t dt + A0 &1(t) = -0.000218 e (485) In differential notations, Equation 485 takes the form (D + A0) l(t) = -0.000218 e+t (486) Equation 486 has the following general solution: -A t t 0 e A C1 1(487) l(t) = C1 e - 0.0(0218 D + A0' A- (487) where C1 is an arbitrary constant determined from the requirement that at t = 0, g1 = =0. Again, by using Equation 483, Equation 487 can be written in the following form: (t) = Ce -(+A 0.000218 et et] ( t < 1), A 1 (488) 0.002 18 Ce _+ A 0 The information obtained from Equation 488 about the actual error c(t) depends on the choice of the value of the parameter AO. From the physics of the problem and a consideration of Equation 481, we can conclude that A0 is < 0, since x(t) < 0. Figure 32 shows different solutions for different reasonable values of A. It indicates that the value of the error Eok < 1.7 El[x(1)] = 0.00037, which gives results correct to the first three digits. The method of successive approximation prescribed above is programmed for the 707 digital computer, which solves the problem by the steps shown in Figure 33. The results are given in Table V which also contains the exact solution so that one can calculate the true errors resulting from the approximations. The true errors are also plotted in Figure 32. The maximum value of the error is 0.000083 which is much less than 0.00037; this agrees with the predicted results. Figure 34, which gives the different approximations of the actual solution, shows that the curve representing the fourth iteration approaches the exact solution of the nonlinear system in such a way that the curves representing the remaining iterations, corresponding to n = 5, n = 6, and 118

Institute of Science and Technology The University of Michigan 2.8 2.6 A =-2 2.4 2.2- 1. ~~~~~~~0 1-0.5 2 - 1.8 -A0 =-.6A1.4 1.2 - ~1.2.3,4.5.6.7 o 8.9 1 t FIGURE 32. A PLOT OF THE ACTUAL ERROR PROPAGATION AND THE THEORETICAL ERROR PROPAGATION FOR DIFFERENT VALUES OF A0 n = 7, almost coincide. Obviously, the curve corresponding to n = 0 gives the solution of the linear system: dt +x= et(489) so that t = O, x = 1. 119

Institute of Science and Technology The Oniversity of Michigan Input Values of Startx At N tmax) n=O0 t =At t=t+At ~IL, t + t te2 -2t No x t) e max' Xn(t) = ~ N Yes n=l t = t + A t=At 12~~~~~~~~E0 k=k+l A~t A t - L~t? No _ = E+ [ ] J+ [X n- 1 (kAt)] Yes Xn(t) = e (e- (At) 1/ 2 + E+ [1/2et] [ 21 (t] t t+ At No m ax Yes FIGURE 33, A TYPICAL FLOW CHART FOR THE SOLUTION BY THE METHOD OF S UC C ESSIVE A PPROXIMA TION 120

CCD 0 TABLE V. NUMERICAL RESULTS FROM THE DIGITAL COMPUTER'* CD Answer Errors 3 4 5~~~~~~~~~~~~~~~~~~~~~~ by the Exact due to RungeN 0 1 2 3 4 5 6 ~~~~~~~~~~~~~~Method of Answer Successive Kutt Approximation x -- -e-t Approximation Method - - 0.05 -0.904837 -0.949086 -0.951137 -0.951235 -0.951239 -0.914 -0951240 -.520 -.529 -.001 -.5292 -0,95240 -.95120 -0.51229-0.00011 -.95129:0 0.10 -0.818731 -0.897050 -0.904311 -0.904823 -0.904855 -0.904857 -0.904857 -0.904857 -0943 -.002 -.087 0.15 -0.740818 -0.844979 -0.859211 -0.860615 -0.860728 -0.860736 -0.860737 -0873 -.678 -.002 -.677 0.20 -0.670320 -0.793685 -0.815666 -0.818463 -0.818742 -0.818766 -0.818768 -0.818768 -0873 -.007 -.170 0.25 -0.606531 -0.743761 -0.773598 -0.778232 -0.778786 -0.778841 -0.778845 -0.778846 -0781 -.004 -.780 0.30 -0.548812 -0.695628 -0.732983 -0.739802 -0.740751 -0.740859 -0.740870 -0.740871 -0.740818 -0005 -.781 0.35 -0.496585 -0.649569 -0.693831 -0.703073 -0.704537 -0.704725 -0.704745 -0.704747 -0.704688 -0005 -0748 0.40 -0.449329 -0.60576~ -0.656162 -0.667958 -0.670048 -0.670345 -0.670381 -0.670385 -0.670320 -0006 -.732 0.45 -0.406570 -0.564308 -0.620003 -0.634384 -0.637195 -0.637634 -0.637691 -0.637698 -0.637628 -0007 -.662 0.50 -0.367879 -0.525239 -0.585373 -0.602292 -0.605897 -0.606508 -0.606594 -0.606605 -0.606531 -0007 -.653 0.55 -0.332871 -0.488545 -0.552285 -0.571625 -0.576076 -0.576888 -0.577011 -0.577028 -0.576950 -0007 -.564 0.60 -0.301194 -0.454178 -0.520740 -0.542337 -0.547661 -0.548701 -0.548869 -0.548892 -0.548812 -0.000080 -0581 0.65 -0.272532 -0.422069 -0.490729 -0.514383 -0.520588 -0.521875 -0.522096 -0.522128 -0.522046 -0.000082 -0524 0.70 -0.246597 -0.392128 -0.462233 -0.487721 -0.494794 -0.496345 -0.496625 -0.496668 -0.496585 -0.000083 -0.496585 0.75 -0.223130 -0.364254 -0.435223 -0.462311 -0.470223 -0.472048 -0.472394 -0.472450 -0.472367 -0.000083 -0426 0.80 -0.201897 -0.338341 -0.409663 -0.438112 -0.446819 -0.448924 -0.449341 -0.449411 -0.449329 -0.000082 -0492 0.85 -0.182684 -0.314277 -0.385509 -0.415085 -0.424533 -0.426916 -0.427409 -0.427495 -0.4,27415 -0.000080 -0.427414 CD 0.90 -0.165299 -0.291951 -0.362712 -0.393190 -0.403316 -0.405973 -0.406543 -0.406647 -0.406570 -0.000077 -0.406569 C 0.95 -0.149569 -0.271253 -0.341220 -0.372384 -0.383120 -0.386043 -0.386692 -0.386815 -0.386741 -0.000074 -0.386741 1.00 -0.135335 -0.2 52075 -0.320978 -0.352629 -0.363903 -0.367077 -0.36780 7 -0.367949 -0.367879 -0.000070 -0.367870 CD U,~ 0.-h -(0 K) 0~~~~~~~~~~~~~~~~~~~~~~~

Institute of Science and Technology The University of Michigan 1.0-.9.8.7-.6-.5 - 4 n4 n 3.3 -n2 n= 1 n 0. I X p 1. ts ~ _I I I I 5,,, I.1.2.3.4.5.6.7.8.9 1.0 FIGURE 34. DIFFERENT APPROXIMATION FOR THE ACTUAL SOLUTION 5.6. COMPARISON WITH SIMILAR METHODS Similar methods of analysis have been developed bIy Stout [36] and Wolf, Ku, and Dietz [42]. These methods require the conversion of a nonlinear differential equation to a nonlinear integral equation of the convolution type. In the notations used thus far, it has the form x(t) = x (t) +{ W(t - u) N(x, 1,. du (490) in which x 0(t) contains the response of the linear portion of the system to both the forcing function and the initial conditions. A step-by-step method of analysis is then used to approximate the integrals. This method results in algebraic equations which, in general, are nonlinear and which must be solved for each discrete value of x k. Stout solved these nonlinear algebraic equations graphically; however, this method does not provide the means for eliminating extraneous solutions which can occur for certain nonlinearities. Wolf, Ku, and Dietz solved these equations analytically and suggested a method for eliminating the extraneous solutions by diff erentiating

Institute of Science and Technology The University of Michigan in which each step in the solution results in a new equation with a new set of initial conditions; this permits the computer to use more accurate rules of integration. Stout was unable to calculate the discretization errors in the solutions because the graphical procedures he used restricted him to working with a nonlinearity which was a function of either x or k, only. The method of Wolf, Ku, and Dietz requires a starting procedure to calculate the first (M - 1) discrete values of the solution necessary to carry on the computations. The numerical analysis method presented herein differs from these other methods by using successive approximations of integral equations of the superposition type, e.g., t x(t) = xo(t) + W(t, u) N(x,,... )du (491) 0 and by including, in the analysis, systems with time-varying parameters. The integrals of different iterations are evaluated by a numerical procedure carried out on a general purpose digital computer. While investigating electric circuits, L. A. Pipes [32] adopted the same method for the approximation of integrals of the convolution type represented by Equation 490. This requires analytic expressions for the nonlinear function so that the integrand increases in complexity even under the best conditions. Thus, the method of Pipes is only suitable for very small nonlinearities which require one or two iterations, and does not provide the means to judge the accuracy of the solution obtained. In this report there is no complexity of integrals since numerical techniques are employed, and the errors committed during the solution can be evaluated. Thus, the soluti on can be obtained to any degree of accuracy. The need for a starting procedure and the possibility of obtaining extraneous solutions are eliminated. Further, the M - point integration rule, described by Wolf, Ku, and Dietz, can also be employed to make possible the use of more accurate rules of integration. 5.7. CONCLUSIONS AND REMARKS A systematic method has been presented for extending the general partition method to the numerical analysis of a broad class of physical system described by higher order differential equations. The forcing function, the nonlinear function, and the time-varying parameters can be specified graphically, tabularly, or by explicit mathematical functions. A means of estimating the errors in the approximate solution has been devised so that the accuracy of the results can be determined and the appropriate interval selected. Although developed for the analysis of

Institute of Science and Technology The University of Michigan An important property of the analysis presented is that, for a given nonlinear term N, the complexity of the analysis does not increase greatly as the order of the system n gets higher. If the system is of a higher order, then the numerical solution of a linear equation of higher order should be found, a relatively simple procedure. This makes the method convenient and practical for the analysis of higher order systems. 6 CONCLUSIONS, CRITICISMS, AND SUGGESTIONS FOR FUTURE STUDY 6.1. SUMMARY The use of the partitioning technique in nonlinear analysis was briefly reviewed, and it was stated that under a broad class of conditions, properly imposed on the linear, nonlinear, and forcing function terms, a class of systems exists for which exact and unique solutions can be found. Mathematically, the class of systems of interest can be represented by an equation, such as Equation 1, which satisfies the conditions of Equations 2 through 4. The types of nonlinearities allowed include algebraic nonlinearities in the form of polynomials of x and its admissible derivatives, and transcendental nonlinearities restricted to functions of x(t) whose inverse has a derivative which is a rational function. In the authors opinion, the present method will prove most useful for nonlinear synthesis rather than nonlinear analysis. A complete synthesis theory could be developed by using the system shown in Figure 6 as a basis. Then, by appropriate partitioning, almost any desired configuration could be obtained for single and multiloop feedback systems.- Arguments both for and against the method of analysis were given in detail. Note, however, that the analysis given in Section 1 is confined to physical systems but does not include all physical systems. This is evident at once from the use of the Lipschitz condition (Equation 4) which is a sufficient condition for uniqueness only. Thus, any other condition which would assure uniqueness could have been selected, and from that an analysis could have been made. Hence, in Section 2, the analysis was extended to study the behavior of a class of systems that can be described by a nonlinear differential equation, with time-varying parameters, of the form of Equation 96. By placing certain restrictions on the linear, nonlinear, and forcing function terms-, the equation was found to have a unique solution which is the limit of a sequence of iterates {x} each of which together with the system so40 lution belong to an Lspace.- Under the conditions given, uniform convergence of the sequence {x}I to the solution x(t) is guaranteed in the mean sense and the ordinary sense. 41The solution. 40

Institute of Science and Technology The University of Michigan x(t) was shown to remain less than some upper-bound function 13(t) which also belongs to L 2 space in the period of interest 0 < t < T. The upper-bound function 1(t), which can easily be obtained, was shown to be a linear combination of some functions in L2[0, T] which are defined from the system restrictions and which depend on the system parameters. Thus, the control engineer can gain further insight into his system by examining the conditions imposed on the system equation; hence, this technique should help solve problems of nonlinear synthesis. A study of the behavior of the admissible system with respect to initial conditions and the system parameters also was presented. The study showed that the system response x(t) is. uniformly continuous with respect to both initial conditions and system parameters. Equation 170 gives an upper bound for the error when the actual solution x(t) is approximated by the n-th iterate x in the sense of mean square. The number of iterations n required for this approxin mation was shown to depend on how the system equation is partitioned. Since the L space is complete, the solution x(t) also will belong to L2 space. This is not true for the case studied by A. A. Wolf [47]. Wolf demonstrated that the limiting process of a sequence of functions of zero order and exponential type can give rise to functions of higher order that are not of exponential type. In Section 3 the partitioning technique and the state-variable approach were used to study the bahavior of the system Equation 96 when it was governed by certain conditions. Under these conditions, the state variables (v1, v2.., vn), which represent the state of the system in the phase space, were shown to belong to L2 space. By studying the definition of asymptotic stability in the sense of Lyapunov, it was found that the system, as restricted, satisfied that definition. Hence, the restrictions placed on the system were sufficient not only to prove uniqueness of a solution that belongs to L2 space, but also to guarantee, in the sense of Lyapunov, asymptotic stability. By selecting two different definitions of the norm, the behavior of the system trajectory was examined for two cases: (1) at any time t and (2) on the average during the interval of interest. In many applications, studying the system trajectory on the average is usually unsatisfactory, since the average value may not exceed a certain specified value whereas the actual system response may assume the undesirable shape given in Figure 23. For this reason an expression for an upper-bound state vector inside which the system will always remain during operation was obtained. The upper,-bound state vector was found to depend only on the restrictions placed on the system; hence, it can be changed to suit the specific application by changing the system parameters. Although matrices are valuable tools for simplifying the

Institute of Science and Technology The University of Michigan le to the form of equations previously analyzed; therefore it was more useful to impose the restrictions on the system equations as they were. (This procedure becomes clearer in the numerical analysis of Section 4.) In Section 4, the analysis presented in Sections 2 and 3 was extended to include systems whose dynamic behavior can be represented by a set of n simultaneous nonlinear equations. Sometimes, it is difficult to reduce these equations to the canonic form of a single-degree-offreedom nonlinear equation of the form studied in Sections 2 and'3. Section 4 showed that a class of multiloop systems exists for which the exact analytic solution can readily be obtained to any degree of accuracy within a prescribed error. This class of systems is defined by the restrictions given in Section 4, and the stability and boundedness of the solution are guaranteed under the assumed conditions that are given there. As an incidental result of the work in that section, the relations between n terminal-pair circuit configurations and multiloop feedback system configurations were investigated and explained. The outputs of the feedback systems were shown to be equivalent to the current response of the circuits, and an example was given to illustrate the ideas and methods of analysis presented. In Section 5, a systematic procedure for extending the general partitioning method to numerical analysis was presented. The forcing functions, the nonlinear functions, and the timevarying parameters can be specified graphically, tabularly, or by explicit mathematical functions. Equations representing the propagation of the discretization error equations in the approximate solutions were developed. This development made it possible to judge the accuracy of the results and to select the appropriate interval. The method is applicable to the analysis of systems described by a single equation, and can be extended to the analysis of systems described by simultaneous equations. A detailed discussion comparing the given method of analysis with similar methods by Stout and Wolf, Ku, and Dietz was given, and arguments were presented in support of the various approaches. 6.2. QUESTIONS FOR FUTURE WORK This report has discussed only a few of the many aspects of the analysis of nonlinear systemns. Let us cite some of the aspects which still should be investigated. (1) Exploration of other criteria to include a liew class of nonlinear systems. For example, the system restrictions can be chosen so that the solution belongs to L space for p > 1. The

Institute of Science and Technology The University of Michigan associated with these norms. They constitute an important class of constraints for the input or the output of a given system, a fact of paramount importance in many engineering problems today. We observe here that when p = 1, p = 2, and p = co, the norm of x(t) is equal, respectively, to the area of x(t), the square root of the energy of x(t), and the maximum magnitude of x(t) in the interval [0, T]. (Actually, when p = o, ix pi is equal to the essential supremum of x(t) in P the interval [0, T].) In order to help the reader follow the analysis when p = 1, it is convenient to make the following assumptions, suggested by H. T. Davis in Reference 11. (a) The function x0(t) is integrable and bounded in the interval 0 < t < T, where x0(t) has the same meaning as in our analysis. (b) The weighted nonlinear function F(t, u, v) satisfies the following two conditions: (i) F(t, u, v) is integrable and bounded, that is, F(t, u, v) < M in a given domain, when M is a finite number. (ii) F(t, u, v) satisfies the Lipschitz condition F(t, u, v F(t, u, v )I < LIv 2 1 2 within its domain of definition. (2) Expansion of the solution x(t), since x(t) E L2[O, T], by a linear combination: cy(t) + c2(p2(t) +... of an orthonormal system of functions qo.(t) over [0, T], and development of methods for finding recurrent relations for the c's and the general solution of these recurrent relations. (3) Development of procedures for solving systems such as those discussed in this report, but with random excitations. (4) Development of a stability theory along other lines of reasoning than those presented here. (5) Extension of the discretiza'tion error equations of Section 5 to include systems described by simultaneous differential equations. (6) Application of partitioning methods to solve systems whose linear parts consist of both lumped and distributed parameters which can be described by a lumped-distributed parameter nonlinear system as follows: L (D, t)x(t) + L (D), t)x(t - u) + N(x) = g(t). 2~~~~~~~~~~~2

Institute of Science and Technology T h e U n iversity of Michigan Appendix A THE EXACT SOLUTION OF EQUATION 30 Given: dx + x + ax = (492) so that t = 0, x = 0. Separating the variables in Equation 492 gives dx 2d- -dt (493) 2 6 ax +6x - Integrating both sides of Equation 493 gives V2+46 1g(8 + 2ax) - 482 + 4a94 log =-t + C 44 V62 + 4a6 (6 + 2ax) + 62 + 2a3 where C is a constant of integration. Applying the initial conditions gives 2 g -V82 + 4ap C=- log - 2 (495) 62+4e3 La 6 + 2aJ Substituting for C in Equation 494 and solving for x gives x= ~~ +a/3_ r 61 +- -d 2+ 2a: at _= 6+.l 6 + 162 + 4a3 _____a -at(4 6 wh _ere where a = + 1/2 V2. +4ac(

Institute of Science and Technology The University of Michigan Appendix B A SUMMARY OF THE L2 SPACE Definition (1) In referring to an L2 space, sometimes called the Hilbert space, we mean the totality of quadratic integrable functions f(t) defined in the closed interval [0, T] for which the norm, or the length I1f fl as given by fff If2 f0 f2(x) dx exists and is finite. Definition (2) The kernel g(t, u) defined in the square (O < t < T, 0 < u < T) belongs to an space if the norm, as given by If g l2 = JTI g2(t, u) dt du exists and is less than a certain constant N Fubines Lemma For the kernel g(t, u) as defined in Definition (2), the functions A(t) (t, u) du and T0 B (u) = g(t, u) dt exist almost everywhere for (0 < t < T, 0 < u < T) and belong to an 2 2 2 space. Also, g = A (v)dv T B (v)dv. ~0 Statement Iff is in L2 and g is in L2, then (f if) are in L2. In general, any finite linear combination of functions that are in L2 belong to L2. Also, If +gJ _< 1ff 11 +g Igf; that is, the triangle inequality holds. Statement If f is in L2 and g is in L2, then the product fg is absolutely integrable and Schwarz's inequality is given byL (fg) du < f2 (u) du g2(u) du Statement 2 If f(x) is integrable over a finite interval, then f(x) also is integrable over a finite interval; but the converse is not true. (Note that this complete statement is not true for infinite intervals.) Theorem Given a sequence {fn} of L2 functions, if tfI} converges in the mean to the L function f, T 2 T 2 that is, nLimn I [fn(t) - f(t)] dt = 0, it always follows that Lim J0(f f n.o0 m___oo 0

Institute of Science and Technology The University of Michigan H. Weyl's Lemma Given a sequence of L2 functions in an interval (0 < t < T), if Lim f = sp(t) (convern~~~~~~~~ 2 n —co n gence in the mean) and Lim f = L(t) (convergence in the ordinary sense), then so(t) = 4(t) n-o n almost everywhere. Appendix C42 THE IMPULSE RESPONSE OF TIME -VARYING LINEAR SYSTEMS AND THE EXTENSION OF THE ANALYSIS OF SECTION 2 TO SITUATIONS IN WHICH THE FORCING FUNCTION IS gl(t) = K(D, t)g(t) WHERE K(D t) IS A LINEAR TIME -VARYING OPERATOR The behavior of linear time-varying systems can be described in several ways. Usually, the behavior is described implicitly by means of a differential equation relating the input g(t) to the output x(t). This equation is, in general, of the form n~~~~~ [an(t)Dn + + a+ 1(t)D + a0(t)]x(t) = [bm(t)Dm +. + bI(t)Dt + b 0tg(t) (498) where D = d/dt and the a's and b's are known functions of time. Equation 498 can be written more compactly as L(D, t)x(t) = K(D, t)g(t) (499) where L(D, t) and K(D t) represent the left-hand and right-hand operators, respectively, of Equation 498. A commonly used form of explicit description of the behavior of system 498 is based on the concept of the impulse response of the system. The impulse response W(t, u) is the solution of the equation: L(D, t)w(t, u) = K(D, t)6(t - u) (500) subject to the homogeneous boundary conditions at t = u or, equivalently, to the following conditions: W(/X)(t~ u) = 0 for X = 0, 1,2,..,(n - 2), and W (nl ) (t, u) = at I#0 where 6(t - u) an (u) represents a unit impulse at t = u. The salient property of W(t, u) is that the response to any given input g(t) can be expressed in terms of W(t, u) and g(t) through the superposition integral 00 x (t) = W (t, u) g(u) du (501) 130~~~~~~~~~~~~~~0

Institute of Science and Technology The University of Michigan If it is assumed that g(t) = 0 with t < 0 and that a physically realizable system, that is W(t, u) is always zero for t < u, Equation 501 reduces to x(t) = W(t, u)g(u) du (502) f: In the analysis presented in Section 2, a special case of Equation 499 was considered, namely, equations of the form L(D, t)x(t) = g(t) (503) By definition, the impulse response of Equation 503 is identically equal to the Green's function of Equation 499, which is denoted by G(t, u). Obviously, G(t, u) should have the same homogeneous boundary conditions given to Equation 500. Note that W(t, u) and G(t, u) are closely related to each other (as shown in the following equations). Treating K(D, t)6(t - u) as the input to the system Equation 503 give( l 00 W(t, u) = G(t, T) [K(D, 7) 6(T -u)d (504) 0co where T is the variable of integration and D = d/dT. Making use of the identity f (T) D(k) 6 (T -u)dT = 1)k f(k) () (505) permits Equation 504 to be -reduced to the desired relation, namely:, W(t, u) = b0(u) G(t, u) - b- 1b(u) [G(t, u)]3 +..+ (-1 dm [b (u) G(t, u)] (506) 0 du I ~~~d m m u so that relation 506 can be expressed more compactly as W(t, u) = K*(D, u) [G(t, u)] (507) where D should be interpreted as d/du, and the operator K*(D, u) is the adjoint of K(D, u). Equation 507 shows immediately that the analysis in Section 2 easily can be extended to a new type of forcing functions of the form: g1(t) = K(D), t) g(t) (508)

Institute of Science and Technology The University of Michigan Notice, also, that Equation 507 leads to another conclusion which is characteristic of timevarying linear systems; namely, that the impulse response W(t, u) of the system Equation 498 can be obtained by making use of the impulse response G(t, u) of the simpler system Equation 503. Appendix D A PROOF FOR EQUATION 133 It was proved that (x2 - x1)2 < c2 a2(t) (509) Substituting Equation 509 into Equation 132 when n = 1 gives 2 22 2 (x - x2) < C a(t) a(u) du (510) 0 Repeating the same procedure using Equation 510 and n = 2 gives 2 2 2Ct Cu2 2 2 A2 (x4 - x3) < C a (t) a(u)du a(z)dz < C2a 2(t) 2 (511) 0 0 And it is known that a (t) du < a (t) du < A and, du < du <A2n J-a (u) du {a (z) dz a2(0) d(0)... n times = n' (512) Then, repeating the last two steps for n = 3, 5,..., gives [Xn+2(t) - Xn+1 (t)] < C a (t) (n) (513) n_2 n( (n)! From Equation 513 it can be concluded that rX n+(t) - xn(t) dt < 2 A2n(514), ~~~~~~~~~~~~(n'.) hence, the result. 132

Institute of Science and Technology The University of Michigan Appendix E UNIQUENESS AND SUCCESSIVE APPROXIMATIONS [, p. 53] In order to show that the convergence of successive approximations does not imply uniqueness, let us consider the following familiar example: dx 1/3 d - x (t -=0o, x =0) (515) t ]1/3 X +1 = x0 + IXn(t) du (516) Therefore, x0 =x1 =x2 =x3 =...n. =0 so that the successive approximations are all zero; hence, they converge to the identically zero solution. On the other hand, the function so(t), defined by 2t1/3 S(t)= (518)3 =3 is another solution which exists to the right of the origin. Also, uniqueness does not imply convergence of successive approximations, as shown by the following example: d x9t [t =0 x =0] (519) where ro (t =O,~ -o< x< +0) 2t (0< t< 1 -o x <O0) 4x 2 ~~~~~~~~~~~~(520) 2 2t (0<t <1It < x< +oo) in the region [0 < t < 1, lx, < oc]; the function f (x, t.) is continuous and is bounded by the constant 2. Therefore, uniqueness is guaranteed since Equation 519 is linear. The successive approximations given by

Institute of Science and Technology The Uhiversity of Michigan become x =0, X2n1(t)=t x (t) =t2 (n= 12 (522) Therefore, the sequence {x (t)} has two cluster values for t 0; hence, the successive approximations do not converge. Also -note that neither of the two convergence subsequences x21 and {X2n} converges to a solution because f ~~~~2 X2' = 2t 0 f(t, t) 2n- 1 (523) x2n = -2t f (t, -t) Appendix F EXTENSION OF THE ANALYSIS TO MORE COMPLEX NONLINEARITIES In the analyses given in Sections 2, 3, and 4, the weighted nonlinear function F was assumed (in its domain of definition) to be a function of t, u, and one independent variable, say x. In order to extend the analysis to include more complex situations, that is, to include the situation in which F assumes the new form (where (in) denotes the m-th derivative of x with respect to the independent variable t, m < n - 1 and n is the order of the system), one needs more complex conditions in order to prove convergence, uniqueness, and other properties. For example, consider the conditions required for an analysis when F contains x and x. In other words, F may assume the form F = Frt, u, x, x ].Thus, we need (1) x 0( (t) F, L 2[O; T], where i = 0 or 1, as. above, (i) denotes the i-th derivative of x 0(t) with respect to the time t, and x0(t) is as defined in previous analysis. (2) In a given domain D [ Ix(t)I < d (t), x 1t S d (t), 0 < u < t< T], the function Flt, u, x(u), x (u), and its first derivative with res'pect to t, may satisfy (a) F (i)[t, u, z19 zi(1)1 - FWi [t u, z 2, K2 K(t, Iu) z 2(k) - k k=0 /It /

Institute of Science and Technology The University of Michigan (3) We also need the constraint conditions: Ix (i) (t) < di(t) qi(t) < $ di(t), where qi(t) is given by qi(t) = XI(i)(t) + C a(t) 2( )Am m=O where..,1/2 C = f n2(t) dt] t ~1/2 a(t) =LJt K2(t, u)duj -1/2 A =fJ a2(t)d(t)],,(i_ (i)l. [t(i)rJ. (u) (1),,1. It is obvious from these conditions that the problem is essentially the same as that which was analyzed in detail in Section 4. Similar conditions are given by Equations 353 through 360. The results, method of analysis, and conclusions are the same. Appendix G A PROOF OF EQUATION 289 To prove that Li x(k) _ (k) 54 Lmx0 M x vk+1(54 where x k)and x~k are given by m

Institute of Science and Technology The Uiversity of Michigan xm(k) =X 0(k)(t) - F(k)[t, u, xm l(u)]du (526) requires two steps: (1) One exists. (2) One should prove tthat x() exists.he m-aD m limit in step (1) satisfies Equation 524. For step (1), one can form the infinite series x(k) + (x2(k) _ x1(k)) + (x3(k) x2(k)) + Xl 2 - Xl + 57 1 2 1 ~3 - 2 to an infinity of terms. Note that x (k) is the m-th partial sum of the above series; thus, by m making use of Equation 288, one gets |xmf' (k) + ak(t) Co + C k(t) (528) k=2 Now, applying the M-test on series 527 shows that as m - co, x (k)has a limit which proves P I ~~~~~~~m step (1). For step (2), subtract Equation 525 from Equation 526 and square the result. This gives [x(k) - vk+2 =( {F(k)[t, u, xm 1(u)] - F(k)[t, Lx-k - vk F21 Kkt, u) Km (u) F t~)Id, uj (529 Thus: F(k) v 2 2 T x(u)12 du Lm - k+U < a~k(t) Ix x1(u) - 2(530) It was proved in Section 2 that Lim x (u) = X(t) (531) Therefore, from Equation 526, it can be concluded that Lim x (k) =v(532) m -0 k+1 wa poedinScton2whrei ws h wn hta ed oifntv()tnst

Institute of Science and Technology The University of Michigan Appendix H A CONTINUOUS FUNCTION WHICH BELONGS TO L 0, AND DOES NOT HAVE A LIMIT AT t = +' A continuous function, which is square integrable on [0, +orc), need not have a zero limit at t = +oo. This can be seen in Figure 35. From Figure 35 we have f(t)dt = 1/4 n=1 In this example the limit does not exist. In order for the limit to exist, the following lemma [50, p. 86] is needed: if f(t) is uniformly continuous4 on [0, +cc), and if f(t) L for p < c p then Lim f(t) = 0; a similar statement holds for a function of more than one variable. t-cc,f(t) 1-. AA FIGURE 35. A CONTINUOUS FUNCTION f(t) WHICH BELONGS TO L2 [0, oo] AND DOES NOT HAVE A LIMIT AT t = +oo

Institute of Science and Technology The University of Michigan Appendix I CONVERGENCE IN THE MEAN AND CONVERGENCE IN THE ORDINARY SENSE The following example shows that convergence almost everywhere does not imply convergence in the ordinary sense [30, p. 169]. Let a sequence {fn(t)} be defined in [0, 1] by the requirements f (t) = n for 0< t <and n ~~~~~~~~n -n = 0, otherwise. For arbitrary t E [0, 1], Lim fn(t) = 0 almost everywhere; but, at the same time, n -co n 1 2 1/n 2 f (t) dt n dt = n. Therefore, f converges to zero almost everywhere, but does not fl2n2 n converge to zero in the mean. The following example shows that convergence in the mean does not imply convergence everywhere [30, p. 96]: 11 ~~~~~~~(k) ()f()t For the half-open interval [0, 1], define the k functions f1( (t)f(t), for every natural number k, according to the relation, 1 forte t i:, Ef.(k) (t) k i ot'k~i Lootherwise in particular, f1 ~(t) = 1 on [0, 1). If all these functions are numbered with indexes 1, 2, 3,,the sequence can be obtained. Now, it can easily be shown that convergence in the mean is zero and that, at the same time, Lim yo (t) = 0 is fulf illed f or no point of t in the interval [0, 1). In f act, if n-ob0 n t0E [0, 1), we can find an i for every k so that and fi k) (t0)1 Inohrwrs0omte o a eg ntesqec fnmes 1t) 2t) 3t)

Institute of Science and Technology The University of Michigan. Appendix J INTERPOLATION, APPROXIMATE DIFFERENTIATION, AND APPROXIMATE INTEGRATIONIn calculating the derivative or the definite integral of a function by means of a set of given values of that function, the function is usually represented by an interpolation formula which can be differentiated or integrated as desired. Briefly, the general roblem of interpolation is to represent a function f(t), known or unknown, in a form po(t), chosen in advance, with the aid of given values which the function f(t) takes for discrete values of the independent variable t. For the purpose of our work the function Sp(t) usually is chosen as a polynomial and the interpolation is known as polynomial interpolation. The justification for replacing a given function by a polynomial rests on a theorem given by Weierstrass in 1885. This theorem states: Every function which is continuous in an interval (a, b) can be represented in that interval to any degree of accuracy by a polynomial p(t). It is possible to find this polynomial p(t) so that If(t) - p(t) I < 8 for every value of t in the interval (a, b) where t is any preassigned positive quantity. This is clarified geometrically in Figure 36. The interpolation formula usually is given in terms of the first, second, and higher order differences which are obtained from a diagonal or horizontal difference table, shown in Tables VI and VII. Ole ~~~~~y f(t).0.0.0 ~~~~~~~y =f(t)~~P(t) _____ _____ _____ t Axis t~~~~a ~~~t b

Institute of Science and Technology The University of Michigan TABLE VI. DIAGONAL DIFFERENCE TABLE tx A2 3 4 5 6 7 8 t X AX A X A X. Ax Ax x x Ax Ax x 0 0 to ~1 A x0 A Ax0 t x2 A2 x0 Ax A 1x t x 2 A2x A x 2 X1 5X0 AX 2 A Ax2 A x A 2 4 6 t 3 X 3 AX2 A XI X 0 3 5 7 t x A Ax A A Axx0 43 232 2 2 6 t4 x Ax A4 A 6 XA 8 Ax4 A x Ax A x1 t5 x5 m 24 m ~3 m62 Ax5 4X 4 2 5x t X AX A X 6 6 5 4 6 2 A 2 t7 x7 A X6 Ax7 t8 x8 From Table VI, the diagonal difference table, we obtain First differences: AxO =x1 - x0, AxI=x2 -X I.. A. A XnI =x -x 0 012 n-i n n-i 2 Second differences: A x = Ax - Ax = x- 2x + x 20 1 0 2 1 0' Ax =Ax2- AxI =x3 -2x2 +x1, 3x2 Third differences: A 30 =A 2X -A 2x0 =x3 -3x + 3x - x 30 2 2 0 3 2 1 0 A1 XIx 2 X x 1 =x4 3x3 +3x2 x1 From Table VII, the horizontal difference table, we obtain First differences: Ax1X = X1I -x1,AIX2- =-XI,1..A IX Xn - Xn Second differences: A2X2 =AIX2 -Ax XI X22x I+ x0 Ax % Ax -rAx =x -3x +3x. -x c

Institute of Science and Technology The University of Michigan TABLE VII. HORIZONTAL DIFFERENCE TABLE [t Ax A x A2x A3x A X AxX 6x A7x A8x ~~ X 1 2 3 4 5 6 7 8 t x 0 0 t x Ax t x Alx A x 2 2 1 2 2 2 t x3 Alx A2x A3x 33 1 3 2 3 3 3 t4 x4 Ax4 Ax4 Ax 4 Ax t5 X5 2lX5 A2 3X5 43X5 A4X 5 5X5 t x Ax Ax Ax Ax Ax Ax 6 6 16 2X6 3X6 4X6 A5X6 A6X6 t7 x7 A2 lX7 3 27 A4 3X7 5 4X7 6 5X7 A 7X7 t x8 A x8 x AX8 A x A X8 A x A x A X8 8 8 1 8 2 8 3 8 4 8 5 8 6 8 7 8 8 8 3 Note that A xI = A3x4 and, in general, the relations between the A's with exponents and those with subscripts are Axk = A xk+ (going forward from xk) xk = kA mk+r or (533) A X = Amxn-m (going backward from xn) where m denotes the order of differences, and k and n denote the number of tabulated values. In the solution of differential equations, the discrete values of the functions for previous times are always known, but those for future times are always unknown. The horizontal difference table shows all orders of differences in terms of the last known values of the function; these differences can be used to compute the next value of the function. For this- reason, the horizontal difference table is more convenient than the diagonal one. In situations where it is necessary to interpolate data at the beginning of the table or near the middle of the tabulated set, one usually uses formulas containing differences given in the diagonal difference table. Thus, for example, Newton's formula for forward interpolation is given by x0 +ux0 +u(u -1) A2x,nh- - 4-,, - A Nu Xr L V A r 2 0

Institute of Science and Technology The University of Michigan where x0 = the value of the function at t = to h = interval width = x1 - x and is assumed uniform, i.e., a constant independent s+l s of x,ands= 0, 1,...,n-1 t - to - u h or t = to + hu h~~~ Anx0 is as explained in Table VI 0 Equation 534 is used when the interpolation is near the point (t, x0 Note that the inter polation function so(t) of Equation 534 is an n-th degree polynomial representation of the original function which is given in terms of the values xO, x1,..., and which is the value of the function for the values t0, t1,., of the independent variable. If interpolation is required near the center of the given tabular set, a central difference formula is usually adopted. For example, we can use Bessel's central difference formula: 2 2 A x A x xx (u ) - 1 x0 (u- 1/2)(u)(u- 1) X X +~~~~uu Ax++- 1) 2 x0 +u x0 +2 2 3 4 4 A4x + A4x u(u - 1)(u + 1)(u - 2) 2 -1 (u - 1/2)(u)(u - 1)(u + 1)(u 2) 5 e+ ALx 4! 2 5+ 2 u (u -)(u + 1) (u - 2)(u + 2)(u - 3) -3 -2 + 6! 2 2 n 2 n u (u -)(u + 1) (u - 2) (u + 2)... (u -n) (u + n - 1) -n + -n+1(5) + ~~~~~~~(2n)! X 2 where the difference appearing in Equation 535 can be obtained by a special consideration of the quantities lying as near as possible to the horizontal line drawn halfway between x 0and x I These quantities are shown in Table VIII. If u = 1/2 in Equation 535, the formula of interpolating to halves can be obtained: x0 +x 1 Ax + Ax A x + Ax A x +A x ____ -1 0 - 2 -1 -3 -2 x 2 8 2 ~ 128 2 1028 2

Institute of Science and Technology The University of Michigan TABLE VIII. CENTRAL DIFFERENCE TABLE 2 3 4 5 6 7 8 X Ax Ax AX Ax Ax Ax A x A x xh x_4 Ax_4 2 x Ax4 X3 -4 3Ax A x -3 X-4 2 4 x AxA X -2 -3 -4 Ax 2 A4 63 X-_ x 2 Ax3 Ax4 3 5x 74 Ax AAA1 Ax2Ax A x-2 1Ax -3 AX- 4 Ax2 2 1 6 0 3 l -2 X,, 5-3 X- 7 x -2A A Ax A~~~~~~~~~ _ 2 4 6 8 / x1 x2 3'" 3 5 7 2 ~~~~~~~~~~4 Ax4 =Ax+1 A 4x + u 1 x u 1 A x 2 2 3x Axn + x+)u2(~)A 2 x2 Ax x$ x3 A~ 3x2 x4 Ax3 Ax4 x5 below (where horizontal difference tables are used): ~0,(t) =;o(tn + hu) = V/(u)

Institute of Science and Technology The University of Michigan The problem of calculating the derivatives of a function by means of a set of given values of that function is solved by representing the function by an interpolation formula such as those given by Equations 534, 535, and 536 and differentiating this formula as many times as desired. The considerations governing the choice of a formula employing differences are the same as those governing interpolation. That is, if the derivative of a function at a point near the beginning of a set of tabular values is desired, use is made of Newton's formula given in Equation 534; but, if the derivative at a point near the end of a table is desired, use is made of Newton's formula given in Equation 537. For points near the middle of the table, a central difference formula should be used, for example, Bessel's formula given in Equation 535. When the above t - to formulas are used to obtain derivatives, one should make the following substitution: u h dx dx du 1 /dx_ andt du dt h du For example, using Newton's formula 534 to calculate the first few derivatives at t = t0 yields 21 /dx 1 [ 2u - 1 2 3u -6u +2 3 dt =Lh Axo + 2! - x0 + 3 x( I- A xo+ A I +...8 /dt x 2 2 3! 6 tdt2 ~ ) +2[0 0 3 Ax At t =t u =0, Equations 538 reduce to the following: 27 Ax x dt Ax0- 2! 3!x (d ) =.[2 x0A3x +..1 (539) 0 Evienlywecanfidx th deiaie3 neatytesm a b ifrnitn h te ie interolatin forulas 144~~~~~~

Institute of Science and TechnologyThe University of Michigan The problem of numerical integration, like that of numerical differentiation, is solved by representing the integrand by an interpolation formula and then integrating the formula between the desired limits. In this way a quadrature formula is derived; the process is known a mechanical quadrature if it is concerned with functions of a single variable, as is the present case. Integrating Newton's formula 534 with the limits of integration taken as t t0and t = t0 + nh or u = 0 and u = n gives: t +nh 0 x(t) dt u (u- 1) 2Xo u(u- 1)(u It h 0 A fl 0 + 3ZA x + du 0 fI or (540) 2 3 2 3 2 A 4 A Ihnxn n n 0 n 3 2 l~ 2 ~+3 2 2: +-~,- n +n +.. From Equation 540, several well-known special formulas can be obtained: the trapezoidal rule, Simpson's third rule, Simpson's three-eighths rule, Weddle's rule, and others. In the case of 2 the trapezoidal rule, for example, n = 1 is assumed and therefore all terms containing A x0and higher differences are rejected, i.e., the curve joining the two points is assumed to be a straight line. Thus, it can be shown that t0+h 0 Similarly, t0+2h= { = [xh x] and so on t0h Adding all such expressions gives integration from t 0 to t0 + nh (when n is even or odd), as f ollow s:

Institute of Science and Technology The University of Michigan where yi = (1, 2, 3,...). For Simpson's one-third rule, let n = 2 and neglect all differences beyond the second in the general formula 540. Then use three ordinates at one time to obtain t0+2h xdt =(x0 + 4x1 +x2) 0 t +4h xdt = (x2 + 4x3+x4), and soon 0o By adding all such expressions, the following formula by Simpson is obtained: t +nh h I=& x(t)dt=[x0 +4(x1 +x3 +. +x n-2) 2(x2 +4 +. 1) +x] to n h-v = Z. (542)ii when0i = 1, 4, 2,..., 4, 2, 1, andeven. when.1, 4, 2,...,y4, 2,1, and neven. If n =6 and all differences beyond 6 are neglected, and the six coordinates (x0, x 1, x) are used at one time, Weddle's rule is obtained. For the first six intervals, this gives 6 t0 ~6h 0 ~~3h {tTo + x, +6x + x + 5x xdt-j(+x x2 x3 4 5 6 For the next six intervals, from x6 to x12 it give s 10(6 x12 12 J' xdt =-1 + 5 7~x8 + 6x 9 + x10 + 5x 11+ 12 6 Finally, adding these expressions gives t0 +n h n x 53 t 0 x(53

Institute of Science and Technology The University of Michigan A computer should have some means of estimating the reliability of every computed result. There are some means of estimating the magnitude of the majority of unavoidable errors. As an example, the error committed in the evaluation of an integral by using the trapezoidal rule, usually known as the truncation error, is given by 2 E <+-x (544) T - 12 m 1 where T = nh = the interval of integration x" (t) = maximum value of the second derivative of the integrand x(t) m A more optimistic estimate of the error is given by the equation 2 ET < + -12 [x'(nh) - x'(O)] (545) 1T 1 where x'(nh) and x'(O) are the first derivatives of the integrand at the end and the beginning of the interval T. Since the first derivative is the integral of the second enh x'(nh) - x'(O) = x"(t) dt = T xav (546) where x" denotes the average value of the second derivative of the integrand in the interval T. av Therefore, Equation 545 can be written 2 E <Th, (547) T <+-f 12 av 2 Equation 545 is much easier to use since it contains only the difference between two easily located first derivatives and not the maximum value of the second derivative, the determination of which is more difficult. Moreover, Equation 545 gives a fairly accurate estimate of the truncation error for reasonable values of T. For these reasons, Equation 545 was used in the calculations. -Of cour se the first derivative can easily be obtained, to any degree of accuracy, from the previous considerations in this appendix. In estimating a given integral by Simpson?'s rule, the truncation errors committed are undoubtedly less than those resulting from use of the trapezoidal rule, but they are more difficult to obtain. It can be shown that the truncation error in Simpson's rule is given by the following equation [35, Chap. 8]: 4 E Th (TV) (48

Institute of Science and Technology The University of Michigan It is not always possible for a computer to have an explicit formula which gives the value of the truncation error committed in the evaluation of a given result. One example may be the Runge-Kutta method. It is true that Bieberback44 has found an expression that provides an upper-bound for the error at a given step in the Runge-Kutta process, but since this estimate depends on quantities that do not appear directly in the computation, some additional separate calculations are required. For this reason, the author suggests the use of the general method of Section 5.2 as a means for calculating truncation errors on the digital computer.

Institute of Science and Technology The University of Michigan REFERENCES 1. Beckenbach, E. F., Modern Mathematics for Engineers, University of California, Engineering Extension Series, McGraw-Hill, New York, N. Y. 1956 p.406. 2. Bilal, A. Y. and Kazda, L. F., "The Analysis of a Certain Class of Nonlinear Time Varying Differential Equations," presented at the Joint Automatic Control Conference at New York University, New York, June 1962, and accepted for publication in the Trans. Am. Inst. Elec. Eng. 3. Birkhoff, G. D. and Kellogg, O. D., "Invariant Points in Function Space," Trans. Am. Math. Soc., 1922, Vol. 23, pp. 96-115. 4. Bochner, S., Lectures on Fourier Integrals, Princeton University Press, Princeton, N. J., 1959, Chap. VIII. 5. Bohn, E. V., "Concerning the Partition Theory and Associated Transform Methods," Proc. I.R.E., July 1961, Vol. 49, No. 7, p. 1215. 6. Bulletin of the National Research Council, Division of Physical Sciences, "Numerical Integration of Differential Equations," No. 92, 1953. 7. Cesari, L., "Problems of Asymptotic Behavior and Stability," Am. Inst. Elec. Eng., Applications and Industry, September 1961, Vol. 80. 8. Coddington, E. A. and Levinson, N. J., Theory of Ordinary Differential Equations, McGraw-Hill, New York, N. Y., 1955. 9. Crandall, N. Y., Engineering Analysis, A Survey of Numerical Procedures, McGraw-Hill, New York, N. Y., 1956. 10. Cunningham, W. J., Introduction to Nonlinear Analysis, McGraw-Hill, New York, N. Y., 1958, p. 122. 11. Davis, H. T., Introduction to Nonlinear Differential and Integral Equations, United States Atomic Energy Commission, 1961. 12. Ingewerson, D. R., "Principal Definitions of Stability," Am. Inst. Elec. Eng. Workshop on Lyapunov's Second Method, September 1960, L. F. Kazda (ed.), The University of Michigan, Ann Arbor, Mich. 13. Ince, E. L., Ordinary Differential Equations, Dover Publications, Inc., New York, 1956, pp. 296-403. 14. Kalman, R. E and Betram, J. E., "Control System Analysis and Design Via the Second Method of Lyapunov," Trans. ASME, Ser. D., June 1960, Vol. 82D, pp. 371-393. 15. Kamke, E., Differentialgleichungen, Lods ungsmethoden und Ldos ungen, Leipzig, 1943, Vol. 1. 16. Kantorovitch, L. V. and Akilov, G. P., Functional Analysis in Normed Spaces, Moscow, 1959. 17. Kaplan, W., Operational Methods for Linear Systems, Addison Wesley, Reading, Mass., 1962, Chap. 8. 18.V- Krnc G. M.(T and Savrachik, T P.E,"An Application oif Fuinctional Analyscis, to an

Institute of Science and Technology The University of Michigan 20. Ku, Y. H., "Analysis of Nonlinear Coupled Circuits II," Trans. Am. Inst. Elec. Eng., September 1955, Vol. 74, pp. 439-43. 21. Ku, Y. H., "Analysis of Nonlinear Coupled Circuits," Trans. Am. Inst. Elec. Eng., 1954, Vol. 73, pp. 626-31. 22. Ku, Y. H. and Wolf, A. A., "Transform-Ensemble Method for the Analysis of Linear and Nonlinear Systems with Random Inputs," Proc. Nat. Electron. Conf., 1960, Vol. 15, pp. 441-55. 23. Ku, Y. H., "Theory of Nonlinear Control," Proc. First Intern. Symp. Autom. Control, Moscow, August 1960. 24. Lalesco, T., Introduction A La Theorie Des Equations Integrals, Paris, A. Herrmann, 1912, p. 17. 25. Lanning, J. H. and Battin, R. W., Random Process in Automatic Control, McGraw-Hill, New York, N. Y., 1956, p. 123. 26. Mikhlin, S. G., Integral Equations and their Applications to Certain Problems in Mathematical Physics and Technology, transl. from the Russian by A. A. Armstrong, Program Press, New York, N. Y., 1957. 27. Milne, W. E., Numerical Solution of Differential Equations, John Wiley and Sons, Inc., New York, 1953. 28. Minorsky, N., Introduction to Nonlinear Mechanics, Edwards Brothers, New York, 1943. 29. Murray, F. J. and Miller, K. S., Existence Theorems for Ordinary Differential Equations, New York University Press, New York, 1959. 30. Natanson, I. P., Theory of Functions of Real Variables, transil. from the Russian by L. F. Boron and E. Hewitt, Fredrickingon Publishing Company, New York, Vol. 1. 31. Olmsted, J. M. H., Real Variables, An Introduction to Theory of Functions, Appleton-Century-Crofts, Inc., New York, 1956, p. 72. 32. Pipes, L. A., "Applications of Integral Equations to the Solution of Nonlinear Electric Circuit Problems," Trans. Am. Inst. Elec. Eng., 1953,Vol. 72, pp. 445-50 33. Poincar'e, H., Les Methodes Nouvelles de la Mechanique Ce'leste, Paris, 1892, 1893, 1899, Vols. 1-3. 34. Richardson, C. H., Introduction to the Calculus of Finite Differences, D. Van Nostrand Company, Inc., New York, 1954. 35. Scarborogh, J. B., Numerical Mathematical Analysis, 4th Edition, Johns Hopkins Press, Baltimore, Md., 1958. 36. Stout, T. M., "A Step-by-Step Method for Transient Analysis of Feedback Systems with One Nonlinear Element," Am. Inst. Elec. Eng. Applications and Industry, January 1957, Vol. 28. 37. Titchmarch, E. C., Theory of Functions, 2nd ed., Oxford University Press, London, 1939, p. 248. 38. Tricomi, F. G., Integral Equations, Interscie-ce Pulihes Ne Yrk and

Institute of Science and Technology The University of Michigan 40. Wiener, N., The Fourier Integral and Certain of its Applications, Dover, New York, 1933. 41. Wolf, A. A., "On the Significance of Transients and Steady State Behavior in Nonlinear Systems," Proc. I.R.E., October 1959, Vol. 47, No. 10, p. 1786. 42. Wolf, A. A., Ku, Y. H., and Dietz, J. H., "Systematic Approximation to the Partition Method," Am. Inst. Elec. Eng., Applications and Industry, July 1960, Vol. 1, p. 49. - 43. Wolf, A. A., "Recent Advances in Analysis and Synthesis of Nonlinear Systems," Am. Inst. Elec. Eng., Applications and Industry, November 1961, Vol. 37, p. 289. 44. Wolf, A. A., "Analysis of Transdental Nonlinear Systems," Trans. Am. Inst. Elec. Eng., Vol. 80, in press 1961. 45. Wolf, A. A., "Generalized Recurrence Relations in the Analysis of Physical Nonlinear Systems," Am. Inst. Elec. Eng., Communication and Electronics, September 1961, p. 383. 46. Wolf, A. A., "Recurrence Relations in the Solution of a Certain Class of Nonlinear Systems," Trans. Am. Inst. Elec. Eng., January 1960, Vol. 78, No. 46, pp. 830-34. 47. Wolf, A. A., "A Mathematical Theory of the Analysis of a Class of Nonlinear Systems," Doctoral dissertation, University of Pennsylvania, Philadelphia, Pa., June 1958. 48. Zaanen, A. C., Linear Analysis, North-Holland Publishing Company, Amsterdam, 1956, p. 86. 49. Zadeh, L. A., "Frequency Analysis of Variable Networks," Proc. I.R.E., 1950, pp. 38, 291-99. 50. Zadeh, L. A., "Determination of the Impulse Response of a Variable Network," J. Appl. Phys., July 1950, Vol. 21, No. 7.

Institute of Science and Technology The University of Michigan PROJECT MICHIGAN DISTRIBUTION LIST 5 1 August 1963-Effective Date ~~~~~Copy No.~ Addressee ~Copy No. Addressee I ~~~~Commanding General ~51-52 Director, U. S. Naval Research Laboratory ~U.~S. Army Electronics Command Washington 25, D. C. ~~Fort Monmouth, New Jersey ~ATTN: Code 2027 ATTN: AMSEL-RD 53 Commanding Officer ~~~~~~~~2-3 Commanding General ~U. S. Navy Ordnance Laboratory ~~U.S. Army Electronics Command ~Corona, California ~Fort Monmoouth, New Jersey ATTN: Library ATTN: AM3SEL-CB 54 Commanding Officer & Director ~~~~~~~4-33 Commanding Officer ~U.S. Navy Electronics Laboratory U.S. Army Electronics R & D Laboratory San Diego 52, California ~Fort Monmouth, New Jersey ~ATTN: Library ATTN: SELRA/ADT 55 Commander, U.S. Naval Ordnance Laboratory ~~~~~~~~34 Commanding General ~White Oak U.S. Army Electronics Proving Ground Silver Spring, Maryland Fort Hoachoca, Arizona ATTN: Technical Library ~~~~~~~~~~~~~~~~ATTN: Technical Library 56-80 DDC (TIPCR) ~~35-38 Director, U.S. Army Engineer ~Cameron Station Geodesy Intelligence & Mapping & D Agency Alexandria, Virginia (35) rAt ntet~ligenceDivisioninia81-86 Aeronautical Systems Division (ASD) Wright-Patterson AFB, Ohio (36) ATTN: Research & Analysis Division Wih-atro FOi on (37) ATTN: Photogrammetry Division (81-83) ATTN: ASD (ASRNOO) (38) ATTN: Strategic Systems Division (84) ATTN: ASD (ASAPR-D) ~~~~~~(ENGGM-O~) ~(85-86) ATTN: ASD (ASRNGE-1) 39 Director, U.S. Army Cold Regions Research & Engineering Laboratory 87-89 Commander, Rome Air Development Center ~~~~P. 0. Box 282 ~~~~~~Griffiss AFB, New York Hanover, New Hampshire Gifs FNwYr (87) ATTN: RAALD 40-41 Director, U.S. Army Engineers (88) ATTN: RAWIC Research & Development Laboratory (89) ATTN: RALSS Fort Belvoir, Virginia ATTN: Technical Documents Center 90-94 Central Intelligence Agency 42 Commanding Officer 2430 E. Street, N.W. U.S. Army Research Officee (Durham) Washington 25, 0. C. Hox CM, Duke Station ATTN: OCR Mail Room Durham, North Carolina ATTN: Chief Information Processing Office SI-Sf Scientific & Technical Information Facility P. 0. Box 5700 43 Assistant Commandant Hethesda, Mryland U.S. Army Air Defense School ATTN: NASA Representative Fort Bliss, Texas 97-Si National Aeronautics & Space Administration 44 Commandant Manned Space Craft Center U. S. Army Engineer School Houston 1, Texas Fort Belvoir, Virginia ~~~~~~~ATTN: Chief, Technical Information Division ATTN: ESSY-L 45 Commanding Officer 99 Cornell Aeronautical Laboratory, Incorporated U. S. Army Combat Development Command Wsigo rjcsOfc Intelligence Agency Falls Church, Virginia Fort Holabird, Baltimore 19, Maryland ATTN: Technical Library 46 Commanding Officer, U. S. Army Electronic Research Unit P. 0. Box 205 100 The Hand Corporation Mountain View, California 1700 Main Street ATTN: Electronic Defense Laboratories Santa Monica, California ATTN: Library 47 U.S. Army Research Liaison Office MIT-Lincoln Laboratory 101 Research Analysis Corporation Lexington 73, Massachusetts 6935 Arlington Road Bethesda, Maryland 48-49 Office of Naval Research Washington 14, D. C. Department of the NavyATNCheIfrainadCnolSsmsDvin 17th & Constitution Avenue, N.W.ATN CheIfrainndotolSsmsDvin Washington 25, D. C. (48) ATTN: Code 463 ~~~~~~~102-103 Cornell Aeronautical Laboratory, Incorporated (48) ATTN: Code 463 4455 Genesee Street (49) ATTN: Code 461 ~~~~~~~~~~Buffalo 21, New York 50 The Hydmrograher ATTN: Librarian

Institute of Science and Technology The University of Michia PROJECT MICHIGAN DISTRIBUTION LIST 5 (Continued) Copy No. Addressee Copy No. Addressee 104 U.S. Naval Photographic 107-109 Jet Propulsion Laboratory Interpretation Center California Institute of Techology 4301 Suitland Road 4800 Oak Grove Drive Washington 25, D. C. Pasadena, California 110-111 Battelle Memorial Institute 505 King Avenue Columbus 1, Ohio 105 Coordinated Science Laboratory ATTN: RACIC University of Illinois Urbana, Illinois VIA: Advanced Research Projects Agcy Washington 25, D. (C. ATTN: Librarian VIA: ONR Resident Representative 112 Bureau of Naval Weapon 605 S. Goodwin Avenue Department of the Navy Urbana, Illinois Washington 25, D. C. ATTN: Code RRRE-3 113 Commanding General U. S. Army Materiel Command 106 The Ohio State University Research Foui Washington 25, D. C. 1314 Kinnear Road ATTN: AMCRD-RS-PEColumbus 12, Ohio ATTN: Security Office 114 Commanding'Officer, U. S. Army Liaison Group, Project C VIA: Commander, Wright Air Development Division The University of Michiga Wright-Patterson AFB, Ohio P. OU Box 618 P. Oo Box 618 ATTN: ASRKSE Ann Arbor, Michigan

+ + F AD Div. 15/2 AD Div. 15/2 Inst. of Science and Technology, U. of Mich., Ann Arbor UNCLASSIFIED Inst. of Science and Technology, U. of Mich., Ann Arbor UNCLASSIFIED THE ANALYSIS OF A CERTAIN CLASS OF NONLINEAR I. Bilal, A. Y. THE ANALYSIS OF A CERTAIN CLASS OF NONLINEAR I. Bilal, A. Y. SYSTEMS, by A. Y. Bilal. Sept. 63. 151 p. incl. illu IL U. S. Army Electronics SYSTEMS, by A. Y. Bilal. Sept. 63. 151 p. incl. illu II. U. S. Army Electronics tables, 50 refs. Command/National tables, 50 refs. Command/National (Report No. 2900-392-T) Science Foundation (Report No. 2900-392-T) Science Foundation (Contract DA-36-039 SC-78801/GP-524) III. Contract DA-36-039 (Contract DA-36-039 SC-78801/GP-524) III. Contract DA-36-039 (Project No. 3D5801001) Unclassified report SC-78801/GP-524 (Project No. 3D5801001) Unclassified report SC-78801/GP-524 In nonlinear analysis, the partitioning technique has been IV. Project No. 3D5801001 In nonlinear analysis, the partitioning technique has been IV. Project No. 3D5801001 used to analyze a certain class of nonlinear systems used to analyze a certain class of nonlinear systems whose dynamic behavior can be represented by a non- whose dynamic behavior can be represented by a nonlinear differential equation with time-varying param- linear differential equation with time-varying parameters. When suitable restrictions were placed on the eters. When suitable restrictions were placed on the linear, nonlinear, and forcing function terms, the sys- linear, nonlinear, and forcing function terms, the system equation presented a unique solution which existed tem equation presented a unique solution which existed to the right of the initial state. The system solution was to the right of the initial state. The system solution was given as a limit of a sequence of Picard iterates {Xn} given as a limit of a sequence of Picard iterates {Xn} which are well defined in a given domain and which be- which are well defined in a given domain and which belong to L2 space. A formula was developed which per- Defense long to L2 space. A formula was developed which per- Defense mits determining the number of iterates necessary for Documentation Center mitsdeterminingthe number of iterates necessaryfor Documentation Center (over) UNCLASSIFIED (over) UNCLASSIFIED ~ + + AD Div. 15/2 AD Div. 15/2 UNCLASSIFIEDUNCLASSIFIED Inst. of Science and Technology, U. of Mich., Ann Arbor UNCLASSIE Inst. of Science and Technology, U. of Mich., Ann Arbor THE ANALYSIS OF A CERTAIN CLASS OF NONLINEAR I. Bilal, A. Y. THE ANALYSIS OF A CERTAIN CLASS OF NONLINEAR I. Bilal, A. Y. SYSTEMS, by A. Y. Bilal. Sept. 63. 151 p. incl. illu II. U. S. Army Electronics SYSTEMS, by A. Y. Bilal Sept. 63. 151 p. incl. illu IL U. S. Army Electronics tables, 50 refs. Command/National tables, 50 refs. Command/National (Report No. 2900-392-T) Science Foundation (Report No. 2900-392-T) Science Foundation (Contract DA-36-039 SC-78801/GP-524) IlI. Contract DA-36-039 (Contract DA-36-039 SC-78801/GP-524) Il. Contract DA-36-039 (Project No. 3D5801001) Unclassified report SC-78801/GP-524 (Project No. 3D5801001) Unclassified report SC-78801/GP-524 In nonlinear analysis, the partitioning technique hasbeen IV. Project No. 3D5801001 In nonlinear analysis, the partitioning technique hasbeen IV. Project No. 3D5801001 used to analyze a certain class of nonlinear systems used to analyze a certain class of nonlinear systems whose dynamic behavior can be represented by a non- whose dynamic behavior can be represented by a nonlinear differential equation with time-varying param- linear differential equation with time-varying parameters. When suitable restrictions were placed on the eters. When suitable restrictions were placed on the linear, nonlinear, and forcing function terms, the sys- linear, nonlinear, and forcing function terms, the system equation presented a unique solution which existed tem equation presented a unique solution which existed to the right of the initial state. The system solutionwas to the right of the initial state. The system solutionwas given as a limit of a sequence of Picard iterates {Xn} given as a limit of a sequence of Picard iterates {Xn} which are well defined in a given domain and which be- which are well defined in a given domain and which belong to L2 space. A formula was developed which per- Defense long to L2 space. A formula was developed which per- Defense mits determining the number of iterates necessary for Documentation Center mitsdeterminingthe number of iterates necessary for Documentation Center (over) UNCLASSIFIED (over) UNCLASSIFIED + + +

AD UNCLASSIFIED AD UNCLASSIFIED the approximation of the solution in the mean square DESCRIPTORS the approximation of the solution in the mean square DESCRIPTORS sense. Within the restrictions imposed on the system, sense. Within the restrictions imposed on the system, the system response was found to be uniformly contin- Nonlinear systems the system response was found to be uniformly contin- Nonlinear systems uous with respect to the initial conditions and the sys- Control uous with respect to the initial conditions and the sys- Control tem parameters. Two different definitions of the norms Analysis tem parameters. Two different definitions of the norms Analysis were selected, and the behavior of the system trajecto- Stability were selected, and the behavior of the system trajecto- Stability ry was investigated for two important cases: (1) at any Synthesis ry was investigated for two important cases: (1) at any Synthesis instant of time and (2) on the average during an interval Differential equations instant of time and (2) on the average during an interval Differential equations of interest. The system (under the given conditions) of interest. The system (under the given conditions) was asymptotically stable in the sense of Lyapunov. was asymptotically stable in the sense of Lyapunov. The analysis was extended to include systems whose be- The analysis was extended to include systems whose behavior is governed by simultaneous nonlinear differen- havior is governed by simultaneous nonlinear differential equations. It was found that a class of multiloop tial equations. It was found that a class of multiloop nonlinear systems exists for which it is possible to find nonlinear systems exists for which it is possible to find the exact analytic solution to any degree of accuracy. the exact analytic solution to any degree of accuracy. For situations in which the forcing function, the non- For situations in which the forcing function, the nonlinear functions, and the time-varying parameters can linear functions, and the time-varying parameters can be defined graphically or tabularly, a systematic method be defined graphically or tabularly, a systematic method has been presented to include numerical analysis. has been presented to include numerical analysis. The method lends itself easily to digital computer UNCLASSIFIED The method lends itself easily to digital computer UNCLASSIFIED calculations. calculations. + AD UNCLASSIFIED AD UNCLASSIFIED the approximation of the solution in the mean square DESCRIPTORS the approximation of the solution in the mean square DESCRIPTORS sense. Within the restrictions imposed on the system, sense. Within the restrictions imposed on the system, the system response was found to be uniformly contin- Nonlinear systems the system response was found to be uniformly contin- Nonlinear systems uous with respect to the initial conditions and the sys- Control uous with respect to the initial conditions and the sys- Control tem parameters. Two different definitions of the norms Analysis tem parameters. Two different definitions of the norms Analysis were selected, and the behavior of the system trajecto- Stability were selected, and the behavior of the system trajecto- Stability ry was investigated for two important cases: (1) at any Synthesis ry was investigated for two important cases: (1) at any Synthesis instant of time and (2) on the average during an interval Differential equations instant of time and (2) on the average during an interval Differential equations of interest. The system (under the given conditions) of interest. The system (under the given conditions) was asymptotically stable in the sense of Lyapunov. was asymptotically stable in the sense of Lyapunov. The analysis was extended to include systems whose be- The analysis was extended to include systems whose behavior is governed by simultaneous nonlinear differen- havior is governed by simultaneous nonlinear differential equations. It was found that a class of multiloop tial equations. It was found that a class of multiloop nonlinear systems exists for which it is possible to find nonlinear systems exists for which it is possible to find the exact analytic solution to any degree of accuracy. the exact analytic solution to any degree of accuracy. For situations in which the forcing function, the non- For situations in which the forcing function, the nonlinear functions, and the time-varying parameters can linear functions, and the time-varying parameters can be defined graphically or tabularly, a systematic method be defined graphically or tabularly, a systematic method has been presented to include numerical analysis. has been presented to include numerical analysis. The method lends itself easily to digital computer UNCLASSIFIED The method lends itself easily to digital computer UNCLASSIFIED calculations. calculations.

+ + + AD Div. 15/2 AD Div. 15/2 UNCLASSIFIED UNCLASSIFIED Inst. of Science and Technology, U. of Mich., Ann Arbor Inst. of Science and Technology, U. of Mich., Ann Arbor THE ANALYSIS OF A CERTAIN CLASS OF NONLINEAR I. Bilal, A. Y. THE ANALYSIS OF A CERTAIN CLASS OF NONLINEAR I. Bilal, A. Y. SYSTEMS, by A. Y. Bilal. Sept. 63. 151 p. incl. illu II. U. S. Army Electronics SYSTEMS, by A. Y. Bilal. Sept. 63. 151 p. incl. illu, II. U.S. Army Electronics tables, 50 refs. Command/National tables, 50 refs. Command/National (Report No. 2900-392-T) Science Foundation (Report No. 2900-392-T) Science Foundation (Contract DA-36-039 SC-78801/GP-524) III. Contract DA-36-039 (Contract DA-36-039 SC-78801/GP-524) III. Contract DA-36-039 (Project No. 3D5801001) Unclassified report SC-78801/GP-524 (Project No. 3D5801001) Unclassified report SC-78801/GP-524 In nonlinear analysis, the partitioning technique has been IV. Project No. 3D5801001 In nonlinear analysis, the partitioning technique as been IV. Project No. 3D5801001 used to analyze a certain class of nonlinear systems used to analyze a certain class of nonlinear systems whose dynamic behavior can be represented by a non- whose dynamic behavior can be represented by a nonlinear differential equation with time-varying param- linear differential equation with time-varying parameters. When suitable restrictions were placed on the eters. When suitable restrictions were placed on the linear, nonlinear, and forcing function terms, the sys- linear, nonlinear, and forcing function terms, the system equation presented a unique solution which existed tem equation presented a unique solution which existed to the right of the initial state. The system solutionwas to the right of the initial state. The system solutionwas given as a limit of a sequence of Picard iterates {Xn} given as a limit of a sequence of Picard iterates {Xn} which are well defined in a given domain and which be- which are well defined in a given domain and which belong to L2 space. A formula was developed which per- Defense long to L2 space. A formula was developed which per- Defense mitsdeterminingthe number of iterates necessaryfor Documentation Center mits determining the number of iterates necessary for Documentation Center (over) UNCLASSIFIED (over) UNCLASSIFIED + + + AD Div. 15/2 AD Div. 15/2 UNCLASSIFIED UNCLASSIFIED Inst. of Science and Technology, U. of Mich., Ann Arbor UNCLASSIFInst. of Science and Technology, U. of Mich., Ann Arbor THE ANALYSIS OF A CERTAIN CLASS OF NONLINEAR I. Bilal, A. Y. THE ANALYSIS OF A CERTAIN CLASS OF NONLINEAR I. Bilal, A. Y. SYSTEMS, by A. Y. Bilal. Sept. 63. 151 p. incl. illu, II. U. S. Army Electronics SYSTEMS, by A. Y. Bilal. Sept. 63. 151 p. incl. illu, II. U. S. Army Electronics tables, 50 refs. Command/National tables, 50 refs. Command/National (Report No. 2900-392-T) Science Foundation (Report No. 2900-392-T) Science Foundation (Contract DA-36-039 SC-78801/GP-524) III. Contract DA-36-039 (Contract DA-36-039 SC-78801/GP-524) III. Contract DA-36-039 (Project No. 3D5801001) Unclassified report SC-78801/GP-524 (Project No. 3D5801001) Unclassified report SC-78801/GP-524 In nonlinear analysis, the partitioning technique has been IV. Project No. 3D5801001 In nonlinear nalysis, the partitioning technique hasbeen IV. Project No. 3D5801001 used to analyze a certain class of nonlinear systems used to analyze a certain class of nonlinear systems whose dynamic behavior can be represented by a non- whose dynamic behavior can be represented by a nonlinear differential equation with time-varying param- linear differential equation with time-varying parameters. When suitable restrictions were placed on the eters. When suitable restrictions were placed on the linear, nonlinear, and forcing function terms, the sys- linear, nonlinear, and forcing function terms, the system equation presented a unique solution which existed tem equation presented a unique solution which existed to the right of the initial state. The system solutionwas to the right of the initial state. The system solutionwas given as a limit of a sequence of Picard iterates {Xn} given as a limit of a sequence of Picard iterates {Xn} which are well defined in a given domain and which be- which are well defined in a given domain and which belong to L2 space. A formula was developed which per- Defense long to L2 space. A formula was developed which per- Defense mits determining the number of iterates necessary for Documentation Center mits determining the number of iterates necessary for Documentation Center (over) UNCLASSIFIED (over) UNCLASSIFIED + + +

AD UNCLASSIFIED AD UNCLASSIFIED the approximation of the solution in the mean square DESCRIPTORS the approximation of the solution in the mean square DESCRIPTORS sense. Within the restrictions imposed on the system, sense. Within the restrictions imposed on the system, the system response was found to be uniformly contin- Nonlinear systems the system response was found to be uniformly contin- Nonlinear systems uous with respect to the initial conditions and the sys- Control uous with respect to the initial conditions and the sys-Control tem parameters. Two different definitions of the norms Analysis tem parameters. Two different definitions of the norms Analysis were selected, and the behavior of the system trajecto- Stability were selected, and the behavior of the system trajecto- Stability ry was investigated for two important cases: (1) at any Synthesis ry was investigated for two important cases: (1) at any Synthesis instant of time and (2) on the average during an interval Differential equations instant of time and (2) on the average during an interval Differential equations of interest. The system (under the given conditions) of interest. The system (under the given conditions) was asymptotically stable in the sense of Lyapunov. was asymptotically stable in the sense of Lyapunov. The analysis was extended to include systems whose be- The analysis was extended to include systems whose behavior is governed by simultaneous nonlinear differen- havior is governed by simultaneous nonlinear differential equations. It was found that a class of multiloop tial equations. It was found that a class of multiloop nonlinear systems exists for which it is possible to find nonlinear systems exists for which it is possible to find the exact analytic solution to any degree of accuracy. the exact analytic solution to any degree of accuracy. C) For situations in which the forcing function, the non- For situations in which the forcing function, the non- ( linear functions, and the time-varying parameters can linear functions, and the time-varying parameters can O be defined graphically or tabularly, a systematic method be defined graphically or tabularly, a systematic method -- has been presented to include numerical analysis. has been presented to include numerical analysis. ( 1 The method lends itself easily to digital computer UNCLASSIFIED The method lends itself easily to digital computer UNCLASSFIED o calculations. calculations. AD UNCLASSIFIED AD UNCLASSIFIED the approximation of the solution in the mean square DESCRIPTORS the approximation of the solution in the mean square DESCRIPTORS sense. Within the restrictions imposed on the system, sense. Within the restrictions imposed on the system, CA) the system response was found to be uniformly contin- Nonlinear systems the system response was found to be uniformly contin- Nonlinear systems uous with respect to the initial conditions and the sys- Control uous with respect to the initial conditions and the sys-Control tem parameters. Two different definitions of the norms Analysis tem parameters. Two different definitions of the norms Analysis were selected, and the behavior of the system trajecto- Stability were selected, and the behavior of the system trajecto- Stability ry was investigated for two important cases: (1) at any Synthesis ry was investigated for two important cases: (1) at any Synthesis instant of time and (2) on the average during an interval Differential equations instant of time and (2) on the average during an interval Differential equations of interest. The system (under the given conditions) of interest. The system (under the given conditions) was asymptotically stable in the sense of Lyapunov. was asymptotically stable in the sense of Lyapunov. The analysis was extended to include systems whose be- The analysis was extended to include systems whose behavior is governed by simultaneous nonlinear differen- havior is governed by simultaneous nonlinear differential equations. It was found that a class of multiloop tial equations. It was found that a class of multiloop nonlinear systems exists for which it is possible to find nonlinear systems exists for which it is possible to find the exact analytic solution to any degree of accuracy. the exact analytic solution to any degree of accuracy. For situations in which the forcing function, the non- For situations in which the forcing function, the nonlinear functions, and the time-varying parameters can linear functions, and the time-varying parameters can be defined graphically or tabularly, a systematic method be defined graphically or tabularly, a systematic method has been presented to include numerical analysis. has been presented to include numerical analysis. The method lends itself easily to digital computer UNCLASSIFIED The method lends itself easily to digital computer UNCLASSIFIED calculations. calculations.