2900-91 -T Report of Project MICHIGAN A STOCHASTIC METHOD OF SOLVING PARTIAL DIFFERENTIAL EQUATIONS USING AN ELECTRONIC ANALOG COMPUTER KUEI CHUANG LOUIS F. KAZDA THOMAS WINDEKNECHT June 1960 THE UNIV E RSI T Y OF MICHIGAN Ann Arbor, Michigan

I;, `, Distribution control of Project MICHIGAN Reports has been delegated by the U. S. Army Signal Corps to: Commanding Officer U. S. Army Liaison Group Project MICHIGAN Willow Run Laboratories Ypsilanti, Michigan It is requested that information or inquiry concerning distribution of reports be addressed accordingly. The work reported herein was carried on by Willow Run Laboratories for the U. S. Army Signal Corps under Project MICHIGAN, contract No. DA-36-039 SC-78801, and for Wright Air Development Division under Contract No. AF 33(600)38019. University contract administration is provided to Willow Run Laboratories through The University of Michigan Research Institute.

The University of Michigan Willow Run Laboratories PREFACE Project MICHIGAN is a continuing research and development program for advancing the Army's long-range combat-surveillance and target-acquisition capabilities. The program is carried out by a full-time Willow Run Laboratories staff of specialists in the fields of physics, engineering, mathematics, and psychology, by members of the teaching faculty, by graduate students, and by other research groups and laboratories of The University of Michigan. The emphasis of the Project is upon basic and applied research in radar, infrared, acoustics, seismics, information processing and display, navigation and guidance for aerial platforms, and systems concepts. Particular attention is given to all-weather, long-range, high-resolution sensory and location techniques, and to evaluations of systems and equipments both through simulation and by means of laboratory and field tests. Project MICHIGAN was established at The University of Michigan in 1953. It is sponsored by the U. S. Army Combat Surveillance Agency of the U. S. Army Signal Corps. The Project constitutes a major portion of the diversified program of research conducted by the Willow Run Laboratories 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 Technical Director Project MICHIGAN iii * *

The University of Michigan Willow Run Laboratories CONTENTS Unclassified Preface...................... iii List of Figures.............................vi Abstract................................ 1 1. Introduction............................. 1 2. Stochastic Processes and Related Differential Equations.... 2.1. Markoff Process 2. 2. Discrete Random Walk with Absorbing Boundary 2. 3. Passage to the Differential Equation and Dirichlet Problem in the One-Dimensional Case 2. 4. Petrowsky Theorem 2. 5. Application of the Petrowsky Theorem 3. Two-Dimensional Partial Differential Equations and Stochastic Processes in Electric Circuits................... 3.1. Elliptic and Parabolic Differential Equations 3. 2. The Adjoint Diffusion Equations 3. 3. Electric Circuits and Diffusion Processes 3.4. Two-Dimensional Generalized Dirichlet Problem and Two Independent Circuits 4. Physical Simulation of Diffusion Motion for Problem Solving....... 4.1. Requirements of a Computing System Using the Random-Walk Method 4. 2. Analog-Computer Simulation of the Random Walk 4. 3. Detection of Boundary Absorptions 4.4. Generation of Boundary Values 4. 5. Accumulation of Boundary Values 4. 6. Noise Sources 4. 7. A Complete Dirichlet Problem Solver 4. 8. Solution Time 4. 9. Laboratory Equipment 5. Experimental Results.............. 5.1. Experimental Procedures 5. 2. Laboratory Results: One-Dimensional Case 5. 3. Laboratory Results: Two-Dimensional Case 5. 4. Laboratory Results: Discussion 5. 5. Conclusions 5.6. Extensions Appendix: Boundary-Value Generation with Coordinate Axes Rotated.... 2 2 3 5 6 7 8 8 9 10 13. 15 15 16 19 19 22 22 23 25 26.30 30 31 35 39 41 42.44 References...............48 Bibliography.............................. 48 Distribution List.............................49 V

The University of Michigan Willow Run Laboratories FIGURES Unclassified 1. Analog-Computer Circuit Simulation of a Random Walk.......... 16 2. Power Spectrum of a Physical Noise Source............... 17 3. Analog-Computer Circuit Simulating Equation 61............ 17 4. Graphical Example of a Two-Dimensional Function............ 20 5. Decision Scheme with Function Generators for Generation of O(C)..... 21 6. System for Using the Stochastic Method to Solve the Dirichlet Problem.. 24 7. Experimental Equipment..................... 27 8. Experimental Phototube-Trigger-Reset Circuit............. 29 9. Stochastic Solution of LaPlace Equation, One-Dimensional Case...... 32 10. Determination of'Spectral Density of Noise Generator.......... 34 11. Stochastic Solution of Boundary-Value Problem, One-Dimensional Case.......................35 12. Stochastic Solution of Boundary-Value Problem, One-Dimensional Case, for Two Values of K............. 36 13. Stochastic Solution of LaPlace Equation with Square Boundary....... 37 14. Stochastic Solution of LaPlace Equation with Circular Boundary...... 38 15. Stochastic Solution of Boundary-Value Problem, Equation 82....... 39 16. Stochastic Solution of Boundary-Value Problem, Equation 83......... 40 17. Coordinate Transformation of a Boundary-Value Problem, Two-Dimensional Case......................45 vi

A STOCHASTIC METHOD OF SOLVING PARTIAL DIFFERENTIAL EQUATIONS USING AN ELECTRONIC ANALOG COMPUTER* ABSTRACT This report presents an analog-computer method for solving a general class of two-dimensional boundary-value problems in which only a small number of operational amplifiers is required. The method simulates the boundary-value problem under consideration as a stochastic model by a technique commonly designated the "Monte Carlo method. " In this case, the stochastic method becomes a tool for the numerical analysis of the partial differential equation. Although the Monte Carlo method of digital-computer solution of partial differential equations has been extensively studied, its application to the analog computer is believed to be new. The boundary-value problems for which the stochastic method is applicable belong to a family of generalized Dirichlet problems of the form: 2 2 aD f aaf af af 1 + D 2 K1(x1, x2) l K2(Xl, x2) a x a x ax 2 1 2 f(x1, x2) = 0(x1, x2) on the boundary C where K (x1, x2) and K2(x1, x2) are arbitrary functions of x1 and x2, respectively The boundary C is an arbitrary, finite closed curve-a Jordon curve The boundary-value function 0(x1, x2) is a bounded, single-valued, piecewise continuous function of x1 and x2 D1 and D2 are constants 1 2 The analog computer, using the stochastic method, solves the boundary-value problem at a desired point within the region bounded by C. This solution, which can be shown to be a statistically convergent quantity, appears as the mean value of the summation of discrete boundary values that are randomly generated by the computing system. A number of one- and two-dimensional boundary-value problems of the type described, whose analytical solutions are known, have been experimentally solved. A comparison of the experimental and analytical solutions corroborates the validity of this method. 1 INTRODUCTION The fundamental importance of partial differential equations in our physical universe, and the inherent difficulty of solving them analytically, have resulted in a large number of methods *This research was performed in the Feedback Control Systems Laboratory of the Electrical Engineering Department of The University of Michigan. 1

The University of Michigan Willow Run Laboratories The University of Michigan Willow Run Laboratories for obtaining their solutions. Since the analytic techniques available are limited in number, and when applied often lead to cumbersome results, several important methods have been developed which utilize the capabilities of digital and analog computers. Such machine methods vary in their approach to the problem and have relative advantages and disadvantages as to accuracy, generality, complexity, cost, etc. Analog-computer solution of partial differential equations requires special techniques to overcome the computer's inherent limitation, namely, time, the unique independent machine variable. To solve a partial differential equation deterministically on an analog computer, one must replace some of the equation's partial derivatives by finite difference quotients and thus convert the partial differential equation into a system of simultaneous ordinary differential equations. The major drawback of this method, however, is the number of amplifiers and associated equipments needed. For example, for a ten-station problem (i. e., a partial differential equation reduced to ten simultaneous ordinary differential equations), a minimum number of ten amplifiers must be used. For convenience of adjustment of complex boundary conditions, as many as 30 amplifiers may sometimes be preferred. The purpose of this -report is to present a new approach to the analog-computer solution of partial differential equations. The method is based on the intimate relation that exists between partial differential equations and the random process that arises in the analysis of electric circuits subjected to random excitations. Instead of being used to solve a partial differential equation which governs the probability distribution of the randomly excited circuits, the random process becomes a tool for numerical analysis of the partial differential equation. This probabilistic method of solving partial differential equations can be considered a sophisticated simulation method. It is commonly referred as the "Monte Carlo method." Its advantages are that it requires a smaller number of operational amplifiers than the usual differential-difference method, and that the number of amplifiers and associated equipments does not increase appreciably even if the problem to be solved becomes more complex. 2 STOCHASTIC PROCESSES and RELATED DIFFERENTIAL EQUATIONS 2.1. MARKOFF PROCESS A discrete Markoff process can be interpreted as a process in which the occurrence of an event depends only on the occurrence of the event immediately preceding it. For continuous cases the direct generalization of this statement is not exactly true. In such a case, the continuous process can best be defined by the following fundamental equation known as the ChapmanKolmogoroff equation:1 It is possible to write two equations into a single Stieltjes integral (see Reference 1). 2

The University of Michigan Willow Run Laboratories f(s + t, yl) = f(s, T|x)f(t, yI1t)d (la) integrated over all 4 space where f(t, y |x) dy equals the probability of the vector random process y in the region R at R _ time (t + t0) on the condition that, at time to, the vector random process was at point x. For discrete cases this becomes: P(s + t, ylx)= P(s, |)P(t, y ) (lb) sum over all where P(s, 4 |x) = the conditional probability of finding 4 = 4 at time (to + s) on the condition = x at time t = to The above equations may be described as stating that the passage of a random process y from state x at time t = to to y at time (s + t + t0) must occur via some intermediate state 4 at some intermediate time (s + to). 2.2 DISCRETE RANDOM WALK WITH ABSORBING BOUNDARY For the sake of simplicity, we will restrict our discussion in this section to a onedimensional case; the extension to higher-dimensional cases follows (at least heuristically) immediately. Consider a particle which starts a random walk at a point x lying inside an interval [a, b]. The particle has the probability p(x) of moving to (x + h) and the probability q(x) of moving to (x - h). Let P(t, b x) be the probability that the particle starting at x will reach boundary b at time t = (m At), then P(t, b Ix) satisfies the following Chapmann-Kolmogoroff equation: m m ' P(tm, b x) = p(x)P(tm_, b x + h) + q(x)P(tm_, b x - h) (2) The above equation implies that the procedure for the particle to reach the boundary, b, at t = t = m(At) can be broken down into the following steps: (1) at t = t1 the particle moves to x + h, then reaches the boundary, b, at the next (m - ) At time; (2) at t = tl the particle moves to x - h, then reaches the boundary, b, at the next (m - 1) At time. The boundary conditions for the particle to be absorbed into the boundaries, b, are: P(tm, blx -a) = 0, P(tm, blx >b) = 1 (3) 3

The University of Michigan Willow Run Laboratories The University of Michigan Willow Run Laboratories Similarly, the probability that the particle disappears into boundary,a,before it reaches boundary, b,is given by P(tm, a Ix) = p(x)P(tml, a |x + h) + q(x)P(tm_, alx - h) (4) P(tm, ax a) = 1 P(tm, a|x _b) =0 (5) From these two particle-diffusion probabilities, the following theorem can be proved. Theorem. Given a function U(x, t ) defined by the following difference equation, U(x, tm) = p(x)U(x + h, tM-1) + q(x)U(x - h, tm_1) (6) U(x, t) =A for x < a m - (7) U(x, tm) =B for x > b then U(x, tm) AP(tm, a|x) + BP(tm, b|x) (8) is the solution of Equation 6 with the boundary conditions of Equation 7. Proof. If Equation 8 is substituted in Equation 6, since P(tm, a x) and P(tm, b x) are the solutions of Equations 2 and 4, respectively, it becomes evident that Equation 8 satisfies Equation 6. Equation 8 with the boundary conditions stated in Equations 3 and 5 satisfies the boundary conditions stated in Equation 7. Since P(0, b|x) = 0 for x < b P(t1, bx) = p(x) if x = b - h = 0if x < b- h then P(t1, blx) > P(, b|x) and then by inductive methods it can be proved that P(tm+1, b|x) > P(tm, b|x) Since P(tm, b x) is a probability, its upper bound is unity. From function theory: as tM -o, the monotonic sequence P(t, b x) approaches a limit given by lim P(tm, b x) = P(b x). t — o Similarly, m lim P(tm, a|x) = P(alx) t -oo m 4

The University of Michigan Willow Run Laboratories From the above monotonic properties of these absorbing probabilities it follows immediately that lim U(x, tm) = U(x) lim [AP(tm, ax)] + lim [BP(tm, b|x)] t — oo t -oo t -oo mm m m = AP(a x) + BP(b x) (9) Under these conditions Equation 6 becomes U(x) = P(x)U(x + h) + q(x)U(x - h) (10) U(x) = A for x - a k (11) U(x) = B for x b J Physically, Equation 6 implies that the solution of the boundary-value problem described in Equations 10 and 11 is equal to the mean score of particles that start from x, diffuse into boundary a or b, and make scores equal to A at a, B at b, and zero within the interval [a, b]. 2.3. PASSAGE TO THE DIFFERENTIAL EQUATION AND DIRICHLET PROBLEM IN THE ONE-DIMENSIONAL' CASE If the functions p(x) and q(x) are so chosen that p(x) - q(x) = -oK(x) + O(h ), where K(x) is a function of x only and O(h 2) is the higher-order function of h, i. e., h, h,.... then, since p(x) + q(x)= 1, p() = 1 + hK(x) + (h2 (12) q(x) and Equation 10 can be written in terms of difference quotients: U(x + h) + U(x - h) - 2U(x) + ~p(x) - q(x) [U(x + h) - U(x) (13 FO (13) h --- —- 2 -- L q(x)h iJ h If h is small enough, then d U U(x + h) + U(x - h) - 2U(x) 2 2 dx h and dU U(x + h) - U(x) dx h where Ax = h. By using the condition stated in Equation 12 and passing h-0, Equation 13 formally becomes the following differential equation: 2 d2U dU d +K(x) = (14) dx 5

The University of Michigan Willow Run Laboratories The University of Michigan Willow Run Laboratories The problem of finding the solution of Equation 14 with boundary conditions U(a) = A and U(b) = B is called the Dirichlet problem, or the first boundary-value problem, for a onedimensional case. In Section 2. 2 it is shown that the mean score of a discrete random walk is a solution of the difference Equation 10 with the boundary conditions of Equation 11. If h is chosen small enough and the condition of Equation 12 holds true, then the mean score of this limiting discrete random walk will give an approximate solution to the first boundary-value problem in a one-dimensional case. 2.4. PETROWSKY THEOREM Petrowsky in 1933 proved an important theorem about the generalized Dirichlet problem and stochastic process for a two-dimensional case. This theorem is fundamental to the present research, in the sense that it guarantees the convergence of the stochastic random-walk problem and the generalized Dirichlet problem. Only Petrowsky's result will be stated here; the details of his proof are omitted, since they are readily obtainable from his original paper (Reference 2). The theorem states: If F (A |x, y) denotes the conditional probability function that a moving point starting from x, y will reach the inner domain of a point region A and the following limiting relations hold: -lim f ) E( -x)F (dE x, y)-M(x, y) (15) lim fxy(X) (15 lim fX() JE (I - y)F(dE, y) - M (x, ) (16) X-0 E AXY E y lim fxy() E( - x)2 (dE[ x, y) - M (X, y) (17) X- 0 W E F lim fxy() ( - x)( - y)F(dE Ix, y) - Mxy(x, y) (18) X-*0 xyE X-0 xy E yy liO f y(A) JE (- i) FA(dEI x y) — Myy(X, ) (19) and lim J ( - x)2 F(dE Ix, y) — 0 (20) A-0 fxy() E-K xy fxy(X) (- y)2Fx(dE lxx, y)- 0(21) -- 0 E-K xy 6

The University of Michigan Willow Run Laboratories where f (A) is a positive function of X and f (X) - 0 as X-0, Kx is a circle of arbitrary xy xy xy radius around (x, y), and E is the domain of interest, then, as X- 0, the function U(x, y) defined by the relation U(x, y) = J (Q)P(dQ|c, y) (22) Q approaches the solution of the following generalized Dirichlet problem 2 2 2 aU a U aU au au M + 2M + M + 2M + 2M-= (23) xx a2 xy axay yy 2 x3x yay Yay+x-+ ay where M M - M > 0 and M > 0 and M > 0, and with boundary condition U(x, y) = xx yy xy xx yy * (x, y), where q(x, y) is uniformly bounded and piecewise continuous on the boundary C. The function U(x, y) defined above could be considered as the mean tally of a random walk starting from point x, y, where O(Q) is tallied every time it reaches a point Q on the boundary C. For the one-dimensional case, the Petrowsky theorem is the same as for the two-dimensional case, except that U(x, y) becomes U(x) and the corresponding generalized Dirichlet equation is reduced to an ordinary differential equation. 2.5. APPLICATION OF THE PETROWSKY THEOREM The importance of the Petrowsky theorem is that it guarantees the convergence of the stochastically obtained solution of the generalized- Dirichlet problem to the analytical solution given by Equation 23. For instance, in Section 2. 2 it was stated that the mean tally of the discrete random walk will approach the solution of the generalized Dirichlet problem for the one-dimensional case. By applying the Petrowsky theorem to this special case, it can be shown that, as h - 0, the mean tally defined by Equation 8 does approach the solution of the Dirichlet problem in a one-dimensional case (Reference 3). In the use of Petrowsky's theorem, the stochastic process which approaches the solution of a generalized Dirichlet problem is not unique; in other words, any stochastic process that satisfies the conditions stated in Equations 15 to 21 will approach the solution of the generalized Dirichlet problem, Equation 23. 7

The University of Michigan Willow Run Laboratories The University of Michigan Willow Run Laboratories 3 TWO-DIMENSIONAL PARTIAL DIFFERENTIAL EQUATIONS and STOCHASTIC PROCESSES in ELECTRIC CIRCUITS 3.1. ELLIPTIC AND PARABOLIC DIFFERENTIAL EQUATIONS If t = t + At and t = t are substituted into Equation 6, then the following equation m m-i results: U(x, t + At) = p(x)U(x + h, t) + q(x)U(x - h, t) (24) Assuming the following relations, which were given in Section 2. 3, p(x) + q(x) =1 and 2 (25) p(x) - q(x) = 2K(x) + O(h) Equation 24 can be written in the following difference-quotients form: U(x, t + At) - U(x, t) _ q(x) U(x, At + - U(x, t) q(x)[(x + h, t) + U(x - h, t) - 2U(x, t)] h U(x + h, t) - U(x, t) ( + K(x) At(26) From the relations of Equation 25, it is seen that q(x) is of the form q(x) =: - K(x) — O(h2) (27) 2 If At = h and Equation 27 is substituted into Equation 26, then letting At - 0, Equation 26 formally becomes the following differential equation: aU(x, t) 1 a U(x, t) K(x)a(x, t) (28) ~Tt- ^~2~2 — "^ ^ "ax - (28) at -2 2 ax ax for t > 0, U(x, t) =A for x = a; U(x, t) = B for x = b. Equation 28, which is known as the Kolmogoroff first differential equation, is an important parabolic differential equation in mathematical physics in the areas of diffusion theory and the a vU(x, t) Brownian movement of a particle. U( If t) 0, equation 28 reduces to the form: at 2 a 2+ K(x) = 0 ax2 a = ax or dU dU d + K(x) = 0 (29) dx 8

The University of Michigan Willow Run Laboratories Inspection of Equation 29 reveals it to be the same as Equation 14, which means that the elliptic differential equation (Equation 14) and the parabolic differential equation (Equation 28) have an intimate relation. Roughly speaking, a boundary-value problem with an elliptic differential equation can be considered as corresponding to an asymptotic or a steady solution of a parabolic differential equation with consequent emphasis on terminal positions rather than on the numbers of steps required to reach a given position. 3. 2. THE ADJOINT DIFFUSION EQUATIONS For simplicity, a one-dimensional Markoff continuous stochastic process is again considered here. Since the process is a Markoff process, it implies that its conditional probability density function satisfies the Chapmann-Kolmogoroff equation +00 f(y, t + six, to = 0) =J f(4, six, t0 = O)f(y, t |, s)d4 (30) -00 If s = At, then Equation 30 is of the form +o0 f(y, t + Atx) = f(4, At|x)f(y, tl1, At)d4 (31) -o0 If y and t are fixed and f(y, t | 4, At) is expanded into a Taylor series in 4 around x, then Equation 31 becomes +oo00 f(y, t + At x)= X f(4, At x) {f(y, t x, At) 2 af(y, t x, At) x a+x( At) - x) + a f(y, t x, At)( X)2+.. }d (32) ax t ax since +o0 | f(4, At|x) d =1 (33) -oo and by assuming +o00 lim | (4 - x)f(4, At x)d - a(x) (34) At- 0t -oo 1 +00 lim -t (4 - x)2f(, At Ix) d - 2b(x) (35) At-0 -o0 9

The University of Michigan Willow Run Laboratories and +co-(x +) 2 limrt ( - x) f(4, At x) d - 0 (36) At-0 -oo-(x - 4) Equation 32 reduces to the first Kolmogoroff equation: 2 af a f af af(t, x) = b(x) (t, x) + a(x) a(t, x) (37) at 2 aa ax For fixed x, it can be shown that f(y, t Ix) as a function of y and t also satisfies the second Kolmogoroff or Fokker-Planck equation. In addition, Kolmogoroff's first and second equations are adjointed to each other, the second having the form 2 af a2 a at= [b(y)f(t, y)] -.[a(y)f(t, y) (38) ay It is of interest to note that the conditions stated in Equations 34, 35, and 36 are the same as those stated in Equations 17 and 20 (Petrowsky's theorem), except in this case f y(A) is equal to At. The details of the derivation of the two Kolmogoroff equations for the n-dimensional case can be found in Reference 4. These equations have the following forms: f af n n 2 ': = x1,x Ai(t.,, x',. ' X. n) (39) 9t -1 i 1 x2 n ax. ~n 1 1 1 n ax. ax. 1=11 i=1 j=l I j af n n 2 n a at =1 ayi ay j [BY( Y_ _ * ' ] - ay' [Ai(t' Y1' Y2'.. Yn)f] (40) where =. =a ' (lG-x -x)f d... d4 (41) Bi1j At0O 2At | i i j 1 2 n (41) 1~0" n and 1A. = lin d (4 3. ELECTRIC CIRCUITS AND DIFFUSION PROCESSES 3o 3. ELECTRIC CIRCUITS AND DIFFUSION PROCESSES It has been shown by Wang and Uhlenbeck (Reference 5) that, when an electric circuit is subjected to an exciting input voltage of Gaussian white noise, the response current of the circuit undergoes a diffusion motion or the current response is a Markoff stochastic process. 10

The University of Michigan Willow Run Laboratories Consequently, the conditional probability density function of the current satisfies the ChapmannKolmogoroff equation and Kolmogoroff's two differential equations. If the initial condition, x, of the current is fixed, then the conditional probability density function is a function of the current, y, and t. Therefore, as shown in Section 3. 2, the conditional probability density function f(y, t) satisfies the second Kolmogoroff or Fokker-Planck equation. Given an electric circuit which has an absorbing boundary (that is, whenever the current reaches any value corresponding to the value of the boundary the process is terminated), the conditional probability density function will be a function of the initial value, x, and t. Consequently, the conditional probability density function satisfies the first Kolmogoroff equation, which is the adjoint equation of the Fokker-Planck equation. An example will illustrate this idea. Let an RL circuit, which is subjected to a Gaussian white noise, F(t), be described by the following ordinary differential equation: dy + 3y = F(t) (43) dt where: is a constant and F(t) is a Gaussian random white noise that has the following properties: <F(t)> = 0 (44a) av < F(t)F(t + ) > = 2D6(7) (44b) If the conditional probability density function f(y, t |x) as a function of y and t is of interest, then the conditional probability density function can be obtained as a solution of the following Fokker-Planck equation: af a2 br-4 wi-1^a~y~f] (45) at= 2[b(y)f]-ay [a(y)f] ay where b(y) and a(y) are the conditional moments given by Equations 34 and 35. If in a small time interval, At, the value of the quantity y in Equation 43 changes only a very small amount, then Equation 43 can be approximately written in a difference form: t+At Ay + /y At = | F(t) dt (46) t t+At where $ F(t) dt implies that, although the function y changes very little in time interval At, t F(t) does not necessarily change very little. 11

The University of Michigan Willow Run Laboratories From Equation 46, with the properties of F(t) as given by Equations 44a and 44b, the functions a(y) and b(y) are readily obtained: <aY>rav At a(y) = lim At -{Y -/3y (47a) At-0O and 1 <Ay2> 2D b(y)=- l-im- D (47b) The Fokker-Planck equation for this case is of the form: af a 2f 2 at D + (yf (48) If the current y starts from an initial position x and has a fixed absorbing boundary, then the conditional probability function f(y, t x) = f(x, t) is a function of x and t and accordingly satisfies the first Kolmogoroff equation (Equation 37). Since the first Kolmogoroff equation is the adjoint equation of the Fokker-Planck equation given by Equation 48, then Kolmogoroff's first equation for this case is of the form af = Da2f 3x f(49) at ax2 ax For steady state, a = 0 and Equation 49 becomes an ordinary differential equation: dx D O-x^=0 (50) If the resistance in the RL circuit is nonlinear, then Equation 43 is of the form dy + K(y) = F(t) (51) dt where K(y) is a function of y. The first Kolmogoroff equation then has the form 2 =D 2 -K(x) af (52) at 2 ax ax 12

The University of Michigan Willow Run Laboratories 3.4. TWO-DIMENSIONAL GENERALIZED DIRICHLET PROBLEM AND TWO INDEPENDENT CIRCUITS If two RL circuits with nonlinear elements K1(y1) and K2(y2) are given and the corresponding currents of these two circuits are connected to the x and y inputs of an ordinary oscilloscope, then the conditional probability density function f(y1, Y2, t|x1, x2) satisfies one of the two-dimensional Kolmogoroff equations, depending on whether one treats it as a function of x's or y's. An example of the two-dimensional case is given below. Let two independent RL circuits with nonlinear resistances be given. The equations describing the responses of these circuits subjected to two independent Gaussian white noises are dy1 dt+ K1(Y' Y2) = F(t) (53a) and dy2 dt K2(Y' Y2) = F2(t) (53b) where K1(y1, y2) and K2(Y1, y2) are the functions of y1 and Y2 and F1(t) and F2(t) are the Gaussian random inputs with the following properties: <F1(t)>a = 0 <F (t)> = 0 2 av <F1(t)Fl(t + 7)> = 2D1 (T) (54) <F2(t)F2(t + T)> = 2D26(r) 2 2 av 2 av 2 av OF 1(t) FZ(t)' = <Fi(t)>av <F2(t)>av The conditional probability density of (Y1, y2) is given by the Fokker-Planck equation: at ay Yj [Bijf] a The B.. and A. functions can be calculated using Equations 41 and 42. 1 13

The University of Michigan Willow Run Laboratories For Equations 53a and b and the condition stated by Equation 54, the B.. and A. functions Ij 1 are of the following form: <A 2> 1 lim 1 av 11 2 At-O At 1 <Ay > 1 lim < 1 Y2 av 12 2 At-O At = <A 2> B lim 2 av D (56) 22 2 At-0 At D2 (56) TAy > lirm I av A = At-O-r At = K1(Y Y2) lim <2 >av 2 = At —O At =K2(Yl' Y2) If the conditional probability function, f, is considered as the function of x1, x2 and t, then the f(xl, x2, t) satisfies the adjoint Fokker-Planck or first Kolmogoroff equation. In this instance the following equation results: af a2 f a-f af af D la +D -ax Kl(Y ) Y K (Y' Y) - a (57) a a a 1 t For steady state, Equation 57 reduces to an elliptic partial differential equation. Given a boundary, C, in two-dimensional space and a boundary function O(C), the following procedure is prescribed: The boundary, C, is placed on the face of the oscilloscope whose horizontal and vertical plates are driven by the random voltages y1 and y2. The oscilloscope beam starting from point (x1, x2) within boundary C performs a random walk which eventually reaches boundary C at a point C.. The random walk is then terminated; the boundary point C. is re1 1 corded; the oscilloscope beam is reset to the point (x1, x2). A second random walk is then initiated. After a large number of repetitions of this procedure, the mean value f defined by the equation m - E0(Ci) =f 14

The University of Michigan Willow Run Laboratories approaches in the mean the solution of the following two-dimensional Dirichlet problem: 2 2 a f a 2f af af D + D -K(x 1 - K2(x1, =0 (58) D1 2 + 22 2- ( x x2) ax (58) 1 2122. where f = b(x1, x2) on the boundary C and x1 and x2 are the initial position of the cathode-ray spot. 4 PHYSICAL SIMULATION of DIFFUSION MOTION for PROBLEM SOLVING 4.1. REQUIREMENTS OF A COMPUTING SYSTEM USING THE RANDOM-WALK METHOD Based on the theory presented in Sections 2 and 3, the solution of a Dirichlet problem of the type described by Equation 58 can be obtained through stochastic computing methods, if the following requirements are fulfilled in a physical simulation: (1) Two RL circuits of the form dy dt + K(, y) = F(t) (59a) dy2 dt + K2(1 y2) = F2(t) (59b) must be simulated where K1(y1, y ) and K (y1, y2) are equivalent to coupling resistances having some prescribed characteristics which, in general, will be nonlinear. (2) Two independent voltage sources of white noise, F1(t) and F2(t), with adjustable spectral density D1 and D2, respectively, for driving the two RL circuits must be available. (The simulated two-dimensional random walk resulting from requirements (1) and (2), above, can be observed on a cathode-ray oscillograph if yl and Y2, the response to the driving voltages F1(t) and F2(t), are used to drive the horizontal and vertical inputs, respectively, of a cathode-ray oscilloscope.) (3) The prescribed boundary, C, must be introduced in such a way as to make detectable the absorption of the random movement of the oscilloscope beam. (4) A value of f(Q) must be generated to correspond to any point Q on C at which the absorption occurs. (5) A continuous record of successive generations of O(Q) must be kept or, alternatively, successive values of k(Q) must be continuously recorded. Either method will allow the tallying n of a mean value of. O(Q.). i=1 1 15

The University of Michigan Willow Run Laboratories The University of Michigan Willow Run Laboratories (6) At each occurrence of an absorption, the oscilloscope beam must be reset to the initial point (x1, x2) inside the boundary C, the point at which the solution of Equation 58 is desired and at which each random walk of the oscilloscope beam must begin. (7) This cycle, of operation must be repeated continuously until such time as the value of the mean of the summed boundary values 0(Q) ceases to change; that is to say, the - in at Equation 57 tends to zero or the function, f, reaches its stationary solution. 4.2. ANALOG-COMPUTER SIMULATION OF THE RANDOM WALK An accurate and convenient means of generating the random walk in a real-time simulation includes a d-c analog computer for simulation of the RL circuits having requirement (1) above. Equations 59a and b may be simulated by the operational amplifier circuit shown in Figure 1. Appropriate gain factors in this circuit will, of course, allow D to be varied at will. F(t) -dy y(t) dt K(y) go -K(y) FIGURE 1. ANALOG-COMPUTER CIRCUIT SIMULATION OF A RANDOM WALK An additional requirement is that each computer circuit must be repeatedly reset to a particular initial condition, (x1, x2) (which in general has nonzero value), at random intervals of time. A further reason for the use of an electronic analog computer is its convenience for resetting initial conditions. Requirement (2) of Section 4. 1 forces the use of two independent white-noise sources of voltage; this is, physically, an imposing requirement. This implies that the ideal sources have constant spectral density and Gaussian amplitude probability distribution over an infinite bandwidth of frequencies. A closer examination of the type of noise sources needed will reveal that the actual requirement need not be as rigid upon F(t) as the ideal source implies, in view of the choice of the d-c analog computer used for simulation of the RL circuits. Since the RL circuits simulated on the d-c analog computer are low-pass, narrowband circuits, a physical noise source, F*(t), which has a power spectrum of the type shown in Figure 2, can be used. When F*(t) has a uniform spectral density of bandwidth, w, which is 16

The University of Michigan Willow Run Laboratories four or five times wider than the passband of the RL circuits, it can be considered as a whitenoise source. The particular requirement stated for the bandwidth, 1l, necessary for F*(t) to be considered as a white-noise source, is commonly used in engineering applications requiring a white-noise source. For the general case of a nonlinear element in the circuit of Figure 1, however, the adequacy of this requirement has not yet been analytically proven. However, for a linear case, an analytical approach can be made. Consider the case in which K(y) = C2y (60) When a gain factor, g, is inserted in the forward path of Figure 3a, and if Laplace transformations are performed on the elements of the system, Figure 3b is obtained. It follows that D* (() 4 FIGURE 2. POWER SPECTRUM OF A PHYSICAL NOISE SOURCE (a) Y(s) (b) FIGURE 3. ANALOG-COMPUTER CIRCUIT SIMULATING EQUATION 61. (a) In a time domain; (b) in a complex s domain. 17

The University of Michigan Willow Run Laboratories F*(s) - F(s) - T(s) - + / sg where Y*(s), Y(s), F*(s), F(s) are the Laplace transforms for y*(t), y(t), F*(t), and F(t), respectively. And if D y() is the spectral density of y(t) DF(y) is the spectral density of F(t) (62) DFy*() is the spectral density of y*(t) D F(w) is the spectral density of F*(t) then D D 2 DY () = D () =T(jw) 2 g 2 (63) D D 2 #2 DF* F F c + (C2g) Inspection of Equation 63 reveals that, if g is made small enough, practically speaking the energy at and above cc, in both y(t) and y*(t), can be eliminated, thus yielding Dy* - Dy Since doing this has in no way changed the characteristics of y(t) as a response to a perfect white-noise source, y*(t) should have the characteristics necessary to solve a Dirichlet problem of the type proposed. However, there is a price that must be paid for limiting the gain of the RL circuits. The gain factor, g, which can be defined as the open-loop gain of the circuit of Figure 3, consists of only an integrator. If the open loop is driven by an "error" voltage, E(t), the following results: y(t) =gJE(t) dt = JE(t)[g dt] = JE(T) dr (64) where t A real time and T = gt = machine time. This illustrates that g can equally well be regarded as a time-scale factor for the integrator. It can be seen that, if g is required to be small so that D - Dy, there is an upper bound placed on machine time, which differs from real time by the factor g. This limits the speed with which the random walk will reach the boundary, C, on the average and, hence, the average rate at which the experimental boundary values 0(Q) can be accumulated. The statistical convergence to the solution of any Dirichlet problem by the random-walk method is inversely proportional to the square root of the number of times the beam crosses the boundary, C, which in turn depends on the bandwidth of the noise source; therefore, for a given convergence to the solution, a minimum fixed time is required. Nevertheless, because of the foregoing, it appears that, for some simple functions of K1(yl) or K2(y2), at least a random walk can be physically generated that will have approximately the desired characteristics. 18

The University of Michigan Willow Run Laboratories 4.3. DETECTION OF BOUNDARY ABSORPTIONS The requirement that the absorption of a random walk by the boundary, C, be detectable can be conveniently accomplished by placing on the oscilloscope face an opaque mask which represents the region, R, bounded by absorbing boundary, C. One or more photosensitive circuits, when appropriately shielded from any light except the oscilloscope beam, will signal the exit of the random walk from the region, R. If this technique is used, C can be completely general. The mask of the region, R, will have width and height in the units of the x1 and x2 of Equation 58. Since the oscilloscope beam is driven by computer voltages representing x1 and x2, then by the mere use of the mask for detection, two important scale factors, S1 and S2, will have been implied. S1 relates x1 of Equation 58 to the corresponding x1 voltage, and S2 11 i z relates the x2 of Equation 58 to the x2 voltage. S1 and S2 are factors in the closed-loop gain of their respective RL circuits. Thus, to be correct, they must be included in the gain, g, in Figure 3 and in Equations 61 through 63, in order to estimate the degree to which D -* D and Dy - D. For example, Equation 63, l Y2 2 Y2 amended, becomes Dy*(W) = - - *M DF* (65) y 2 C2) F 4.4. GENERATION OF BOUNDARY VALUES In general, it will be necessary to generate q(Q) at every occurrence of an absorption by C of the random walk, Q being the point of absorption on C. A detailed study into the ways in which a detecting device can be instrumented reveals that the complexity of necessary electronic equipment increases as the complexity of C and 0(C) increases. A detecting device for a completely general C and O(C) would probably be very involved. One possible technique for instrumenting such a device is described in this section; this technique could be appropriately used for generating any arbitrary single-valued O(C), where C satisfies the following conditions: (1) C is a Jordan curve. (2) C can be divided into two parts by a dividing line of the form x2 = mx1 + b such that the arc length (which is usually a function of x1 and x2) of each part is a single-valued function of the length of the dividing line. 19

The University of Michigan Willow Run Laboratories For example, if a line x2 = a can be found which divides C into two parts so that C = C1(xl) for x2 > a and C = C2(x1) for x2 < a, where C1(x1) and C2(x1) are single-valued functions of x1, then the condition is satisfied.2 It is important to point out that the choice of a dividing line which will cause the condition to be satisfied is not, in general, unique. The choice of a dividing line must be left to the ingenuity of the computer operator. However, in the search for an appropriate dividing line, the following property is useful. Of all the line segments of a particular slope, mi, which begin and terminate on C, the only possible dividing line (if it exists) which will cause C to satisfy the condition will be the segment of unique maximum length. In Figure 4, for example, two different boundaries are tested with a dividing line of the form x2 = a to see whether the condition is satisfied. In part a of the figure, the condition is satisfied; in part b it is not satisfied. X 22L Cl(xl) X22 Cl(xl) C2(x,) C2(xl) xX X1 (a) (b) FIGURE 4. GRAPHICAL EXAMPLE OF A TWO-DIMENSIONAL FUNCTION. (a) A singlevalued function. C = Cl(xl), for x2 > a; C = C2(xl), for x2 < a, where C1 and C2 are singlevalued functions of x1. (b) A multi-valued function. C = Cl(xl), for x2 > a; C = C2(x1), for x2 < a, where C1 and C2 are not single-valued functions of x1. If the condition imposed on C is satisfied, then it is possible to approximate the function 0(C) with two diode function generators. In Figure 4, the corresponding boundary-value function, +(C), is of the form: 0(C) = [C1(x1)] for x2 > a q(C) = 0[C2(x2)] for x2 < a The functions 0[C1(x1)] and 0[C2(x2)] can be approximated in straight segments by driving the two diode function generators by a voltage proportional to x1. It should be noted, however, that the boundary-value function is not uniquely defined by the variable x1 alone, but 2The use of a dividing line of more general form than x1 = constant or, x2 = constant somewhat complicates the computational setup. Details are included in the appendix. 20

The University of Michigan Willow Run Laboratories depends also on whether x2 is greater than or less than a. In other words, either O(C) = p [C1(x1)] or 0 (C) = 9 [C2(x2)] is properly the boundary-value function, depending on whether the argument of the decision function, D(x2 - a), is greater than or less than zero. It can be seen that some decision scheme has to be incorporated with the function generators for generation of q(C). One possible scheme is developed below. An opaque divider passing through x2 = a, which is placed perpendicular to the oscilloscope face, separates two phototubes for boundary-absorption detection. When an absorption occurs along any portion of C, a signal is sensed by a phototube, which opens a gate or a related device depending upon the voltage generated from the diode function generator. Thus, at each absorption, a 0(C) is generated by one or the other of the function generators. This operation is illustrated in Figure 5. Random Walk Bounded Region Mask CRO i -- - -, < /!. ~" 'y ----\oyl(t) Phototubes - Y2y(t) Trigger Trigger Pulse Pulse Circuit Circuit 0 [Cl(xl)- Diode Function Gate Generator FU 5 D S D(x -a) _\ _[C2 (XI{)] Diode Function Gate_-_ Generator FIGURE 5. DECISION SCHEME WITH FUNCTION GENERATORS FOR GENERATION OF O(C) 21

The University of Michigan Willow Run Laboratories The University of Michigan Willow Run Laboratories The advantage of this scheme is that, from a practical point of view, the restriction on C is mild, i. e., all common geometric shapes satisfy the condition, and only two single-variable function generators would be required. 4.5. ACCUMULATION OF BOUNDARY VALUES Since a mean value of the experimental boundary values is desired, a technique for continuous summation of successive boundary values and a coexistent tally of the number of boundary absorptions is most convenient. From the standpoint of accuracy, these would be best done in digital form. An example of a digital accumulating technique is included in Section 4. 7 as part of a complete random-walk Dirichlet problem solver. 4.6. NOISE SOURCES A low-frequency noise source of Gaussian amplitude distribution having the power spectrum of Figure 2 is a commercially available item. Such devices generally employ a gas tube subjected to a cross-magnetic field as the source, the power of which is distributed along a band of frequencies in the audio range. This power band is then narrowed by a bandpass filter and chopper demodulated to obtain a power spectrum extending from d-c to somewhat less than 100 cps. These devices have outputs, whose amplitude distribution is closely Gaussian, and have a uniform power spectrum over this frequency range. For the use that is described here, however, an additional requirement must be placed on the noise sources, namely that, over a period of days, weeks, and even months, the spectral densities of each source stay constant. This results from the fact that the machine solutions will be functions of D1 and D2 of Equation 58. Although, in general, D1 and D2 can be varied to values other than the actual densities of the noise-generator outputs by introducing scale factors in the analog computer, the noise-generator spectral densities must be constant, at least for the duration of data-taking on a particular problem. Further, if one is to avoid frequent measurement of these generator spectral densities, a method for long-term constancy must be arranged. A noise generator of the type described ordinarily fails to meet the requirement of longterm stability of spectral density for periods greater than a few hours. One possible way to obtain two noise sources with long-term stability is to record the output of a gas-tube noise source simultaneously on two channels of a dual-channel tape recorder equipped with two playback heads whose spacing along the tape can be varied. If the delay time, T', which for a particular tape speed is proportional to the separation between the playback heads, is made large enough, the two noise sources obtained from these playback heads can be considered to be statistically independent. This results because the cross correlation between the two taped sources at any frequency is a decreasing function of T. 22

The University of Michigan Willow Run Laboratories As has been implied, the requirements placed on the low-frequency response of the tape recorder is severe, but no problem occurs if certain recorders employing frequency modulation are used. As an alternative, however, the output of an audio range source could be directly recorded on a standard dual-channel audio recorder, which again is equipped with variably spaced playback heads. If the playback channels are equipped with identical3 narrow passband filters and demodulators, two independent noise sources are again produced. 4.7. A COMPLETE DIRICHLET PROBLEM SOLVER A block-diagram representation of a system that employs the stochastic method described in this paper, which automatically solves the Dirichlet problem, is shown in Figure 6. The types of Dirichlet problems that this system could solve would be limited to those for which: (1) The feedback elements K1(y1, Y2) and K2(yl, Y2) of the analog-computer RL circuits are functions that will cause the RL circuit to exhibit a low-pass characteristic of cutoff frequency w0, which is much lower than w1, the highest frequency over which the corresponding F(t) may be said to have uniform spectral density. (2) The boundary C satisfies the conditions stated in Section 4. 4. (3) The two boundary-value functions resulting from the condition on C can be reasonably approximated by diode function generators. In Figure 6, an additional refinement has been made, namely, the inclusion of the "hold" circuitry. Its purpose is to hold the RL circuits at the point of boundary absorption, Q., during the millisecond period in which the boundary value is being generated and summed. The generation and summing procedure is accomplished by two identical units which consist of a diode function generator, a phantastron delay circuit, a gate, a high-frequency pulse generator, and a pulse counter. For one of these combinations, at the occurrence of a detection signal from the phototube, a triggering pulse is produced which simultaneously actuates the phantastron delay circuit and the hold circuit. The phantastron delay circuit then produces a pulse of duration proportional to the value of an external voltage, namely, the output voltage of the diode function generator. However, because the hold circuit has also been actuated during the period of pulse generation of the phantastron delay circuit, the diode function generator receives a constant input voltage corresponding to Qi during this period of time. Therefore, the function generator output is a constant voltage proportional to O(Qi), and this causes the phantastron output pulse to have duration proportional to 0(Qi). It will be an advantage if the spectral densities of the two noise sources are the same. Otherwise, the unpleasant task of "matching" the sources would result if a problem were specified such that D1 = D2. 1 2i 23

h) Reset Hold HighFrequency Analog - Diode Pulser At-1 Vz Control Computer Function Circuitry Generator I __3 t-I --- —— l Pulse r"^ g Phantastron- Linear Counter J c~~~~~~~ g * Delay Control Gate x2 Initial g Condition Electronic.2 Monostable Hold Oscillator Cir Circuit Dual-itil Channel, o Condition r 0L T VAnalog- | | | Diod Variable - Dela, Bounded otot C t Computer Pulse Tape Region Disa Reset Counter Noise Source Mask Right Trigger Circuit 3 Ciruitry y1 hototubel Circuit C Feqe Reset Electronic Li _Monostable g Oscillator iru 0 Circuit xF Initial S F N Condition 0 ---PROBLEM ~iw,4 — Phantastron Control Linear Pulse 0__ t[ Delay Gate Cunter Analog- Diode -Computer - FunctionC -High-7 Circuitry Y Generator Frequency t PulserJ Hold Reset FIGURE 6. SYSTEM FOR USING THE STOCHASTIC METHOD TO SOLVE THE DIRICHLET PROBLEM -I ID I3t CD 3 (D -9 I, 0 0 -4 -n O0 I3 Q 0 -9 O CD u <

The University of Michigan Willow Run Laboratories This phantastron pulse opens the gate for a time proportional to O(Qi) and thereby allows a number of high-frequency pulses, again proportional to O(Qi), to reach the pulse counter. The pulse-counter total, therefore, is at all times proportional to: n 1 [ 1( 1)] 2> or [c2(x1 )] X2 <a The total of pulse counter 3 (Figure 6) will at all times be equal to the number of resets that have occurred, which is the same as the number of 0(Qi) values that have accumulated. The mean of the sum of these accumulated boundary values is equal to: A x (total of pulse counter 1) + B x (total of pulse counter 2) (66) total of pulse counter 3 The proportionality constants, A and B, must, of course, be calibrated. The reset of the oscilloscope beam to the point (x1, x2) after each absorption is accomplished by the "reset" circuit. Again, there is in the reset operation a suitable delay for boundary-value generation. 4.8. SOLUTION TIME The basic disadvantage of the system described in Section 4. 7 is the time that is required to obtain the, solution of a problem to a reasonable degree of convergence based on the relatively narrow bandwidths of the modulated gas-tube noise sources previously described, and on the use of a d-c analog computer. A considerable part of an hour would be required to solve a typical problem at one point. The choice of the noise sources and the kind of analog computer described as part of the computing system of Section 4. 7 were limited by the availability of components and the need for economy. However, nothing mentioned so far would prevent the computation rate from being vastly increased if the following modifications in equipment were made: (1) Using a randomly repetitive analog computer for simulation of the RL circuits. (2) Using white-noise sources having wide bandwidths, for example, which can be obtained by tape-recording a gas-tube source which is not demodulated.4 4Some such physical sources exist and have constant spectral density over approximately the same frequency band as a repetitive computer. 25

The University of Michigan Willow Run Laboratories (3) Redesigning the detection and boundary-value generation equipment to require less time for performing their functions. The computing rate, if the above modifications were made, would probably be limited by the operating speed obtainable from the equipment in modification (3). Unfortunately, time has not permitted an adequate investigation of the feasibility of a rapid computing system employing a randomly repetitive computer. 4.9. LABORATORY EQUIPMENT In this section, the experimental equipment used to demonstrate the principles of this stochastic method to solve partial differential equations is discussed. The reader should understand that the purpose of these procedures was to demonstrate the feasibility of the method, not to develop an elaborate computing machine. Features highly desirable in a working computing machine, notably, rapid rate of computation and, to a lesser degree, accuracy of computation have been sacrificed in order that available laboratory equipment might be utilized wherever possible. The computing system used to obtain the results in Section 5 is shown in Figure 7. It is similar to, but is somewhat less elaborate than, the one described in Section 4. 7. The following additional restrictions were applied to all experimental problems considered: [C1(x1)] =T1 =constant x2 > a C2(X1) = T2 = constant x2 < a (67a) or alternatively, [c'(x 2) = T' =constant x > b [C (x2)] = T = constant x < b (67b) These restrictions of course, eliminated the need for function generators, holding circuits, and the equipment described for the boundary-value summing that was presented in Section 4.7. A detection by one phototube of an absorption caused a count on one scalar, and a detection by the other phototube caused a count on the other scalar; therefore, the mean value of the boundary values was merely T1 x (counts of scalar 1) + T2 x (counts of scalar 2) (68) (counts of scalar 1) + (counts of scalar 2) 26

The University of Michigan Willow Run Laboratories Reset (a) From Computer (b) FIGURE 7. EXPERIMENTAL EQUIPMENT. (a) Block diagram; (b) Physical detail. 27

The University of Michigan Willow Run Laboratories The d-c analog computer used was a Sterling Model LM-10, a small, inexpensive, tenamplifier, electronic differential analyzer type. Its amplifiers are unstabilized, but after a brief warmup period exhibit relatively little drift. The reset operation is controlled by a relay in each integrator circuit. The primary noise source was the low-frequency Gaussian noise generator, Pace Model 201 A. This uses a regulated thyratron tube source and employs chopper demodulation to obtain the output. The spectral density of the device is stated by the manufacturer to be constant (~0. 1 db) from 0 to 35 cps and Gaussian distributed in amplitude to ~1%. The spectral 2 density of the device at an average setting was measured to be 0. 96 volt / cps. For the two-dimensional case, the Pace noise-generator output was recorded on an Ampex Model FR 107 tape recorder. This elaborate recorder employs frequency modulation, as previously mentioned, and thereby can reproduce low frequencies extending down to d-c. It was equipped with two sets of seven recording heads and had seven playback heads the position of which was changeable to allow a wide variation of calibrated delays between channels. Only one of four available forward tape speeds was used, namely, 3. 75 ips, which allowed a playing time of more than an hour before rewinding. The delay time used between playback heads was about 2 seconds. The circuit shown in schematic form in Figure 8 was designed and built for detecting boundary absorptions, triggering the two scalars, and resetting the computer. The latter operation was done by driving the computer relays as slaves to the relay shown in Figure 8. While this method gave a slower solution time than others considered, it proved to be more convenient because it required no internal changes in the computer5 and was reasonably dependable. The upper and lower portions of this circuit, when excited by the phototube, generate a single, negative triggering pulse that contains a very steep leading edge. This was necessary because the counting devices used (Tracerlab Model SC-19 utility scalars) had to be used with sensitivities prefixed for nuclear measurements. The remaining portion (reset control) has a recovery-time constant of about 0. 5 second, which was necessary to allow a complete reset of the computer integrators before the next random walk was initiated. Boundary-region masks were fashioned out of cardboard and placed on the oscilloscope (Tektronics Type 532) face, with a plastic piece showing the coordinate axes. The entire system was capable of a modest 2000 random walks per hour with no appreciable error in the detection circuit. 6The computer had to be kept available for other uses during most of the experimental work. 28

CD C +300 (CD ~I I1P37 1 1M 1 Phototube 2 12AX7 22M/ I 1M - _ 0 cq -^ F 150 p/Af *^ 3^.3K 2D21[- '[-II-{ j P ', To Pulse Counter 1 001 — 1 6c I C ^ 4 M< 0.001 1 CD - t +300 r > C1 't O 0 +3000 iL _ _ _ 0 r-f '-_ L - s -loo 16SL7 l. -- 1%~?- 12AU7 2 2 2. 2M 22M im iM 32. 2M - "Reset 1M ^! c~ '^^\^/ ^' ^ ^~i ^ "^ky +300 "? 1 00 1P37 l aflM ^2\ +.loo "Reset" \ 1P37 1 6SL7 cq im LI -100 Relay Phototube 8 EXR2 IMET 12 AX7 - 22M 1 2 -'w 150 ppf '^ i! i- 2D21 - To Pulse Counter 2 0 oo~~~~~ 0-^/$0 =>CI~~~~~~~ C>OC SrT 1_oo -100o FIGURE 8. EXPERIMENTAL PHOTOTUBE-TRIGGER-RESET CIRCUIT O O CD O~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~a

The University of Michigan Willow Run Laboratories The University of Michigan Willow Run Laboratories 5 EXPERIMENTAL RESULTS 5.1. EXPERIMENTAL PROCEDURES The technique used to demonstrate the feasibility of this stochastic method for solving partial differential equations was to solve a number of representative one- and two-dimensional boundary-value problems by the stochastic method, problems in which a graphical comparison of the machine solutions with the known analytical solutions of these problems could be made. In order to compare graphically, the solutions found by the machine and analytical methods, it was necessary to (1) restrict the problems to be considered to those with simple boundary conditions whose analytical solutions can be found, and (2) to solve each problem by the machine method at many points within the bounded region R (for example along a straight line in R). Procedure (2) allows the machine solutions at points along such a straight line to be plotted against the length of the straight line for visual comparison with a continuous curve representing the analytical solution along the same straight line in R. An additional restriction on the kinds of boundary-value problems to be considered experimentally is the restriction on the (C) functions of problems considered. The restrictions stated in Equation 67 were motivated by the simplicity and economy of experimental equipment available to obtain the solutions by the machine method. The stochastic method described herein is probably of more value for solving two- rather than one-dimensional boundary-value problems, because most one-dimensional boundary-value problems, of the type described in Equation 69 below, can be solved in other ways, for example, by analytical methods or by a more conventional use of the analog computer than the one described here. Nevertheless, one-dimensional boundary-value problems can be solved by the method suggested by this report, as has been shown in Section 3, and can be solved with less laboratory equipment than is required for the two-dimensional case; i. e., only one noise generator and one RL circuit are needed for the one-dimensional case. A general statement of the one-dimensional boundary-value problem that can be solved by machine method is as follows: 2f df D 2 -K(x2) dX = 0 dx2 2 2 (69) f(a) = Aa < x2 < b f(b)=B B where a, b, A, and B are constants. 30

The University of Michigan Willow Run Laboratories The differential equation in this case reduces to an ordinary differential equation, and the boundary of the problem consists merely of two points, a and b. The random walk which begins at x2for the one-dimensional case occurs along the straight line connecting a and b (the x2 axis) and can terminate only by absorption at a or b. The simplified scheme and the equipment needed for generating boundary values described in Section 4. 9 for the two-dimensional case can be appropriately used for the one-dimensional case as well, since each boundary value generated is either one or the other of the two numbers A or B. The initial phase of the experimental work of this research was devoted to the comparison of machine and analytical solutions of some representative one-dimensional boundary-value problems of the type described in Equation 69. This was done for the following reasons: (1) A qualitative evaluation of the accuracy of the average machine solution at a point in the one-dimensional case, as a function of the number of random walks used to obtain the solution, was desired. (2) The performance of the experimental phototube-trigger-reset circuit of Figure 8 was tested. (3) The adequacy of the available noise generator was evaluated without complicating this evaluation by tape recording the generator output to obtain two independent noise sources, as would be required for solving two-dimensional problems. Following the attainment of these objectives by comparison of machine and analytical solutions in the one-dimensional case, the additional laboratory equipment needed for solving two-dimensional boundary-value problems was acquired and several two-dimensional problems were treated. Results of these procedures, with graphs (Figures 9-16), are presented below. 5.2. LABORATORY RESULTS: ONE-DIMENSIONAL CASE The first boundary-value problem considered experimentally involves the one-dimensional Laplace equation 2 df =0 (70) dx2 The comparison of the machine solution and its corresponding analytical solution for the particular boundary-value problem solved appears in Figure 9. The number of random walks per experimental point indicated in Figure 9, also indicated on subsequent graphs, refers to the number of boundary values accumulated to find the solution at each experimental point on the graph. 31

The University of Michigan Willow Run Laboratories f -100 d2f 2 dx dx2 ~ -= 0 — 80 f(-1)= 0 f(+l) 100 - - 60 _-/ ANALYTICAL SOLUTION 40 ~ * EXPERIMENTAL DATA -20 -_____ 1 XA< CX~.//Catt I ~lli2/IL -1.0 -0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1.0 X2 FIGURE 9. STOCHASTIC SOLUTION OF LAPLACE EQUATION,ONE-DIMENSIONAL CASE. 700-900 random walks per experimental point. Inspection of Equation 69, reveals the Laplace equation, for which K(x) = 0, as the only case of this problem for which the solution is independent of the value of the input spectral density, D. In order to plot the analytical solutions of boundary-value problems of the type of Equation 69, other than the special case above, the value of the spectral density of the noise generator, Dnoise generator' must be known. This results from the fact that the analytical solutions are generally functions of the ratio K/ D where D, the actual input spectral density of the RL circuit, is related to D. by certain known scale factors of the RL circuit. noise generator The spectral density of the noise generator can be directly measured by several known methods, but, in this instance, an indirect method was used which involved the solution of the following one-dimensional boundary-value problem by the machine method. 2 d f df D Kx2 dx2 f(0) = 0 f(1) =1 0 < x < 1 32

The University of Michigan Willow Run Laboratories which has the analytical solution 1 exp (- K ) f(x 2) = (72) 2 K-ep(-|) If in Equation 72, K >> 1, this solution can be expressed approximately as _K>>1 (73) f(x2)1 - exp (-x2) (73) The procedure for finding the spectral density was first to set up the stochastic computing machine to solve the boundary-value problem described in Equation 71 and then to arrange the value of K to vary in the RL circuit used. Machine solutions of the boundary-value problem of Equation 71 in the range 0. 7 < x2 < 1 were found for various positive values of K until, for a particular positive value of K, the solutions indicated that f(x2) 1 (0.7 <x2 <1) (74) This condition assures that, for this value of K, >> 1, as can be seen from Equation 72. When the condition of Equation 74 was realized, K was fixed and its value recorded. Machine solutions for other values of x2 were then found, so that the point x20, at which f(x20) = 1 - 0. 636 (75) 20 e could be approximately located graphically from a plot of the machine solutions as a function of x2. From a comparison of Equation 73 and 75, it can be seen that, at x20, D = 1 or D=Kx20 (76) The spectral density of the noise generator, Dnoise generator' is related to D by known scale factors in the RL circuit used, and thus can be computed directly. The graphical procedure involved in finding D is shown in Figure 10, from which D 0. 96 volt / cps noise generator - 33

The University of Michigan Willow Run Laboratories I I I d2_f df D 2 + K d d2 2 f(0) = 0 — ]f(l)= i- { - u 1.0 --- 0. f 0. 0. 0. 8 6 4 _ I K -D x2 = 1 1v A -11114 allK4. K D2 1 -e -K/D 1-e K= 2.0 - NOISE GENER EXPERIMENTAL CURVE EXPERIMENTAL DATA EXPERIMENTAL DATA./ / 0 / ATC D I I= 0 9 V I ID = 0.96 VOLT /C 'PS 2 --- -- -- -— {{jJ< — -- -- -- - -- -- -- -- — L - { X __ ___ — 0 0.1 0.*2 0.3 0.4 0.5 x2 0.6 0.7 0.8 0.9 1. 0 FIGURE 10. DETERMINATION OF SPECTRAL DENSITY OF NOISE GENERATOR. 300 random walks per experimental point. The second one-dimensional boundary-value problem to be considered involves the ordinary differential equation, d-f df 2 dx2-K 22 (77) The graphical comparison of machine solutions and analytical solutions of a particular boundary-value problem having this differential equation appears in Figure 11. In this case, the machine and analytical solutions are compared for three different values of K, where the input spectral density, D, has the same value throughout. The analytical solutions that appear in Figure 11 were plotted, using a value for D corresponding to D = 0. 96volt /cps, noise generator which was found by the procedure outlined above. The agreement between the analytical and machine solutions in Figure 11 is good, particularly so because the spectral density of the noise generator was known to change somewhat during the course of data taking. The data taking in this case occurred on several different days, for several hours each day. 34

The University of Michigan Willow Run Laboratories 1. 0. f 0. I I I d2f df D -K =0 2 dx2 f(0)= 0 - - K -0. 5 o _-8 _f(1) =1 _ D = 0. 96 C - __~~ __ __ __ ^~ _ 0. 5 __ _ ~4~~~~-' '' K=:5. 0 2 _' ANALYTICAL SOLUTION - EXPERIMENTAL DATA I II I I Ii 0. 0. 0 0.1 0. 2 0.3 0.4 0.5 X2 0.6 0.7 0.8 0.9 1.0 FIGURE 11. STOCHASTIC SOLUTION OF BOUNDARY-VALUE PROBLEM, ONE-DIMENSIONAL CASE. 300-400 random walks per experimental point. The final boundary-value problem solutions in the one-dimensional case which are included for comparison involve a differential equation of the type d2f df D + Kx-=O 2 dx dx2 2 d22 (78) The machine and analytical solutions for two distinct values of K appear in Figure 12. D, again, has the same value for both sets of solutions. As in the second case, the analytical solutions were plotted with a value of D found by relating the experimental value of Dnoise generator 2 0.96 volt /cps, using scale factors of the RL circuit. 5.3. LABORATORY RESULTS: TWO-DIMENSIONAL CASE The multichannel FM tape recorder was arranged to provide two independent noise sources by recording the output of the gas-tube noise generator simultaneously on two of the recording channels. This was done in such a way as to make the spectral density of the noise on each tape channel approximately the same. Subsequently, a procedure similar to that used to find D, described in Section 5. 2, was used to measure the spectral densities of the noise generator 35

The University of Michigan Willow Run Laboratories e dz.1f 0.60.7 0.8 0.9 1.0 X2 2 ~_/ x2 -4 2 The first two-dimensional boundary-value problem which was solved by the machine method with the tape-recorded noise sources involved the Laplace equation - --- - - - J e dz - 2fa 2f a 2f(79) f = 0 A — AA (79) ax 2 with a square boundary. The comparison of the machine and analytical solutions of the boundary-value problem appears in Figure 13. Comparing the Laplace equation of Equation 79 with the general partial differential equation of the Dirichlet problem reveals that if in Equation 58 D =D (80a) K 1 (1)=S K2(x2)OO 0 (80b) This value is much smeraller than D because the noise-voltage level was o enoise generator lt deliberately reduced during the tape-recording process so as not to exceed the maximum input 36 36

The University of Michigan Willow Run Laboratories f(x1, 0) = 0 x1 1 60 40 - - - ANALYTICAL SOLUTION --- -- -- —._ —* --- EXPERIMENTAL DATA 20 --- 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.C X 2 FIGURE 13. STOCHASTIC SOLUTION OF LAPLACE EQUATION WITH SQUARE BOUNDARY. 1200 random walks per experimental point. then the Laplace equation results. As a result of the restrictions imposed by Equation 80 the solutions of boundary-value problems involving the Laplace equation are independent of the value of D1 = D2. The machine solutions shown in Figure 13 were obtained by using two identical RL circuits on the computer; the circuits were driven by two independent noise sources which had the same spectral density. The omission of machine solutions near the boundary in Figure 13 is explained in Section 5. 4 as part of a more general discussion of experimental results. Figure 14 shows a comparison between the machine and analytical solutions of the Laplace equation with a circular boundary. The Laplace equation for this case is more conveniently stated in terms of the polar coordinates (r, 8) as 2 af a f 1 a 2 V f- + = (81) 2 rar 2 2 ) 37

The University of Michigan Willow Run Laboratories FIGURE 14. STOCHASTIC SOLUTION OF LAPLACE EQUATION WITH CIRCULAR BOUNDARY. 2200 random walks per experimental point. where r Xl2 + x22 r= x2I x2 O tan-1 ( In Figure 15, a comparison is made between the machine and analytical solutions of a boundary-value problem in which the partial differential equation is /2 2f f a f af D + x2 + K = (82) ax 2 ax 2j ax x1 21 and in which the boundary of the problem is a square. The analytical solution of this problem is a function of K/D and was plotted; the value of D was obtained from the value of the spectral density of the tape-recorded noise sources mentioned above. In Figure 16, a comparison is made between the machine and analytical solutions of a boundary-value problem involving the partial differential equation D a X2 +K ax af 0 (83) ax j x + ax 1 a2/ \ 1 2 I 38

The University of Michigan Willow Run Laboratories 1 2 n7r 1 2 n=l ( cs ct i hl s n~r 1a n- - e / (-1)n] (sin n7rx1) L- rX Il - 60 rl 4-M v. a X1 FIGURE 15. STOSCHASTIC SOLUTION OF BOUNDARY-VALUE PROBLEM, EQUATION 82. Two-dimensional case, square boundary. 300 random walks per experimental point. and a square boundary. Again, when the analytical solution was plotted, the value of the spectral density of the tape-recorded noise sources was used for D. When the computing machine was instrumented for solution of the boundary-value problems of Figure 15 and Figure 16, the phototube divider was placed along a diagonal line (x1 = x2) of the square boundary. As mentioned in Footnote 2, the use of a dividing line of a more general character than x1 = constant or x2 = constant requires that the RL circuits used be specifically altered from the usual arrangement in order to solve the desired partial differential equation. The method of properly altering the RL circuits, which is discussed in detail in the appendix, was used to obtain the machine solutions of Figure 15 and Figure 16. These solutions constitute an interesting corroboration of the particular procedure described in the appendix as well as the over-all theory of the stochastic method. 5.4. LABORATORY RESULTS: DISCUSSION Examination of the machine and analytical solutions of various boundary-value problems presented in Section 5. 2 and in Section 5. 3 reveals that, although the machine point solutions 39

The University of Michigan Willow Run Laboratories fxl, (1 - x)]= 4001re2 2 (X 1),n 1 - -- --- - --- -- - - - -. I - - - - - I - - - - - - - - - - 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 x1 X1 FIGURE 16. STOCHASTIC SOLUTION OF BOUNDARY-VALUE PROBLEM, EQUATION 83. Two-dimensional case, square boundary. 300 random walks per experimental point. taken together provide agreement with the analytical solutions, there are point-by-point deviations from the analytical curves. These deviations are caused by the stochastic nature of the solution. In order to discuss the feasibility of the stochastic method as a tool for solving partial differential equations, more must be said of the stochastic convergence of a random process. If the experimental equipment can be said to introduce no error, a machine solution at a point that is the mean value or mathematical expectation of discrete boundary values as they are accumulated by the computing system is a random variable subject to purely statistical variation. Furthermore, the standard deviation of the statistical solution from its true value varies as the inverse square root of the number of boundary values that have accumulated. Therefore, the solutions can be said to be statistically convergent, in the sense that, as the number of accumulated boundary values is increased, the probability is also increased that the deviation of the machine solution from the stationary solution is less than a given amount. The word "convergent" is advisedly used because, for any finite number of accumulated boundary values, there will always be a finite probability, however small, that the deviation of the machine solution will exceed a given amount, i. e., that the solution will not "converge" at all. 40

The University of Michigan Willow Run Laboratories The inherent random deviations of the machine solutions discussed in Sections 5. 2 and 5. 3 are believed to outweigh most errors introduced into the computation results by components of the computing system. There are some exceptions. It is likely that the one-dimensional results, notably those in Figure 11, were influenced by the unstabilized spectral density of the gas-tube noise generator used. In both the one-dimensional and the two-dimensional cases, the solutions at points close to the boundary were erroneously influenced as the result of the limited duty cycle of the reset portion of the phototube-trigger-reset circuit of Figure 8. Because the triggering portions of this circuit have a much shorter recovery time than the reset portion, it was occasionally possible, when the initial conditions placed the oscilloscope beam near a portion of the boundary, for the second of two random walks to be absorbed before the reset portion of the circuit had recovered from the first. When this occurred, the boundary value was erroneously counted a second time. For the machine solution in Figure 13, this difficulty was circumvented by omitting points near the boundary. In later solutions, by reducing the rate of boundary-value accumulation for points near the boundary, the errors caused in this way were reduced sufficiently to allow the machine solutions to be used. 5.5. CONCLUSIONS The question of the feasibility of the stochastic method for solving partial differential equations is largely influenced by the statistical nature of the "convergence" of solutions, as discussed in Section 5. 3. Since this statistical "convergence" is not an absolute or guaranteed convergence as the time of computation is increased, the stochastic method is not a feasible tool for scientific applications requiring a high degree of accuracy. However, the method is potentially useful in engineering applications where: (1) An approximate solution is sufficient. (2) Solution time is not the controlling factor in the choice of a computing method. (3) The bulk, weight, and/or cost of computing equipment must be kept relatively low. The precision7 of a point solution and the time required to obtain the solution are interrelated parameters of the stochastic computing method. For a given computing system, solution time can be reduced only by reducing the precision which can be attributed to the solution; precision can be increased only by increasing the time of computer operation. 7The precision of a point solution can be expressed in terms of a "confidence interval" and a "confidence coefficient." The confidence interval is an interval with the point solution as center within which the stationary value is likely to be situated. The confidence coefficient is a statement of the probability that the stationary value does indeed occur within the confidence interval. 41

The University of Michigan Willow Run Laboratories The feasibility of applying the stochastic method will therefore often depend upon the rate of boundary-value accumulation which can be realized, and can be enhanced by the use of repetitive analog-computer circuitry and the associated equipment mentioned in Section 4. 8. In this regard, it should be emphasized that the precision realized in this research, wherein solutions were obtained by generating only a few hundred boundary values, is encouraging from the standpoint of possible application, and warrants further investigation. A second, more general, observation that deserves some emphasis is this: the analog computer, because of the advantages of real-time simulation, is a natural and convenient tool in the study of random processes. As mentioned in the Introduction, numerical analysis of random processes in the past has been conducted almost exclusively with digital computers. Since the parameters of random processes are statistically varying, the accuracy with which these parameters are controlled and by which they can be measured is not of major importance, particularly if the objective is to study the random processes involved. It therefore seems unnatural to confine the tools of study of random processes to digital computers, wherein the great advantage of accuracy is not fully utilized. 5.6. EXTENSIONS The machine solution of a boundary-value problem can be described as follows: N /T \ f(x1, x2) = n 0(Qn) (84) where T is the number of absorptions occurring at the point Qn on C T is the total number of absorptions occurring during the computation /T The computer data, (, n = 1, 2, 3,..., N), possess a property that is analogous to the Green's function. A set of data, n = 1, 2, 3,..., N), which is found by the machine method for a particular partial differential equation and a particular absorbing boundary, can be used to solve the partial differential equation for any boundary-value function, 0(C). It will merely be necessary to carry out the arithmetic operations indicated in Equation 84 to find the solution for a particular (C); no further use of the stochastic computing machine is required. The Green's function analogy suggests additional applications of the stochastic method. A T set of data, n = 1,2, 3,..., N, can be recorded in some convenient fashion and stored, thus constituting a means of solving a particular set of boundary-value problems of arbitrary k(C), without further need of a stochastic computing machine. There are practical limitations to this procedure. As the absorbing boundary C and function +Q(C) are made more complicated, carrying out the arithmetic operations indicated in Equation 42

The University of Michigan Willow Run Laboratories 84 becomes, in general, unreasonably tedious. At this point, one would have to choose between further use of the stochastic computer method and that of some other computing device. The Green's function property might be well utilized where the boundary-value function, 0(C), is to be varied in order to optimize in some fashion the solution of a boundary-value problem. Again this could be done by manipulating the data, (T) n = 1, 3,..., N), arithmetically after a single run with the stochastic method on the computer. It appears likely that additional investigation, with the present discussion as a base, will lead to the means for solving more general types of partial differential equations than are described here. Among the possibilities is the two-dimensional generalized Poisson's equation (Reference 3). 2 2 a f a f af af D l +D 2 +K -+Ki-^ (85) D1 a2 D2 2 + Kax x = g(xx2) ax ax2 1 1x. 2x 2 1 2 and certain kinds of eigenvalue problems (Reference 6). The underlying theory of the stochastic method discussed in this paper has been extended to boundary-value problems of more than two dimensions (Reference 3). However, the restrictions on the absorbing boundary, C, become considerably greater in this case. There may be some possibilities of extending the machine method to three dimensions. 43

The University of Michigan Willow Run Laboratories Appendix BOUNDARY-VALUE GENERATION with COORDINATE AXES ROTATED The method of boundary-value generation which is described in Section 4. 4 is straightforward if the phototube divider is appropriately placed on a dividing line of the region R which has the form x = constant or x2 = constant (86) However (see Footnote 2), if a dividing line of the form of Equation 86 cannot be used and a more general dividing line of the form x2 = mx1 + b (87) is used, some complication is introduced. The phototube housing for this case must be designed so that the phototubes and divider can be pivoted into any position with respect to the oscilloscope face. In addition, the voltage needed to properly drive the two diode function generators is E = y cos (tan m) + Y2 sin (tan- m) which, for example, may be generated by using a "summing" operational amplifier. During the laboratory work of the present research, an alternative procedure was devised so that the phototube divider (which was fixed in a position parallel to one of the oscillograph coordinate axes) could be used to solve problems for which only a dividing line of the most general form, x2 = mx1 + b, was satisfactory. This procedure is valid only for certain special cases of the general partial differential equation, Equation 58. Consider the effect of placing the boundary mask on the oscilloscope face so that the dividing line, x2 = mx1 + b, is parallel to one of the oscillograph coordinate axes. It is apparent that the boundary-value problem has undergone a coordinate transformation with respect to the oscillograph coordinate axes in which the oscillograph horizontal coordinate axis is referred to as u1 and the vertical axis as u2 (Figure 17). In this figure, the coordinate transformation used is one of pure rotation; that is, the origin of each coordinate system is the same.8 This trans8The transformation can involve translation as well as rotation; however, there is no advantage gained by it for present purposes. 44

The University of Michigan Willow Run Laboratories formation can be described as: x1 = Cu1 - C2u2 2 = C21 + C1U2 (88) where C = cos a 1 C 2 sina 2 U2 X2 X2 { -. U1 = constant / —`R to, x2=mx+ b / I / I ul FIGURE 17. COORDINATE TRANSFORMATION OF A BOUNDARY-VALUE PROBLEM, TWO-DIMENSIONAL CASE 2 2 In addition, since C1 + C2 = 1, the inverse transformation is valid: U1 =Cl + C2X2 2 = -C2X1 + CX2 (89) 45

The University of Michigan Willow Run Laboratories By use of the latter set of coordinate transformation equations, Equation 89, a partial differential equation of the form is transformed, in general, to a form D + D2 af + K (X +K u +K (91) D 1 2 +uXl 1 2f D12'=2C1C2(D2 =a- D KD2t =C D 8+ 2 K = C (x12) + C2K2(x2) K2' =C 1K2(x2) - 02K.jx1) Theis transformed boundary-value problem can be handledal, through use of the computing system of Figure 6 only if Equation 91 reduces to the form of Equation 90 with u1 and u2 substituted for x and x, respectively.f 21 2 2 where D1 = C12D1 + C2 D D 2, 2D 2 ' =C1 + C22D1 K' = K'(u = K'(C + Cx) = C (x) + (92) K11 C1K1(x1) C2K2(x2) =K2(u = C1K2('() + - C1K22K() i(1) Froystem of Figure 6 on of the defquation 91 reduces to the form of Equatitieson 90 with u and the implication ofsubstituted for x1 and x2, respectively. Therefore, the transformation result must be such that D12 = 0 D1 = D2 Kl(x ) = constant (93) LK'X1 K2(x2) = constant LK'x 46

The University of Michigan Willow Run Laboratories where K = constant. If the conditions of Equation 93 are satisfied by the coefficients of the original partial differential equation, two RL circuits of the form dv1 dt + K (Vl) = Gl(t) (94) dv dt + K(v2) = G2(t) where <Gl(t)Gl(t + T)>a = 26D1 ~~~~~~~~~~~~~1 1 av 1"~~~ ~(95) <G2(t)G2(t + T)>av = 26D2 corresponding to Equavtion 91, can be used to solve the original boundary-value problem. corresponding to Equation 91, can be used to solve the original boundary-value problem. 47

The University of Michigan Willow Run Laboratories REFERENCES 1. A. Kolmogoroff, "Uber die Analytischen Methoden in der Wahrscheinlickeitsrechung," Math. Ann., 1931, Vol. 104, pp. 415-458. 2. I. Petrowsky, "Uber des irrfahrtproblem, " Math. Ann., 1934, Vol. 109, pp. 425-444. 3. J. H. Curtiss, "Sampling Methods Applied to Differential and Difference Equations, " Proc. Seminar on Scientific Computation, International Business Machines Corp., New York, N. Y., November 1949. 4. W. Feller, Mathematical Theory of Diffusion, International Congress of Mathematicians, 1950, Vol. II. 5. M. C. Wang and G. E. Uhlenbeck, "On the Theory of the Brownian Motion, II," Rev. Mod. Phys., 1945, Vol. 17, pp. 323-342. 6. M. Donsker and M. Kac, "The Monte Carlo Method and its Applications, " Proc. Seminar on Scientific Computation, International Business Machines Corp., New York, N. Y., December 1949, pp. 74-81. BIBLIOGRAPHY Metropolis, N., and Ulam, S., "The Monte Carlo Method, " J. Am. Statist. Assoc., Sept. 1949, Vol. 44, pp. 335-341. Yowell, E. C., "A Monte Carlo Method of Solving Laplace's Equation," Proc. Seminar on Scientific Computation, International Business Machines Corp., New York, N. Y., December 1949, pp. 89-91. 48

The University of Michigan Willow Run Laboratories DISTRIBUTION LIST 5, PROJECT MICHIGAN REPORTS 1 June 1960-Effective Date Copy No. Addressee 1-2 Commanding General, U. S. Army Combat Surveillance Agency 1124 N. Highland Street Arlington 1, Virginia 3-28 Commanding Officer, U. S. Army Signal Research & Development Laboratories Fort Monmouth; New Jersey ATTN: SIGM/EL-DR 29 Commanding General U. S. Army Electronic Proving Ground Fort Huachuca, Arizona ATTN: SIGPG-DXP 30 Chief of Engineers Department of the Army Washington 25, D. C. ATTN: Research & Development Division 31 Commanding General Quartermaster, Research & Engineering Command, U. S. Army Natick, Massachusetts 32 Chief, Human Factors Research Division, Office of the Chief of Research & Development Department of the Army Washington 25, D. C. 33-34 Commander, Army Rocket & Guided Missile Agency Redstone Arsenal, Alabama ATTN: Technical Library, ORDXR-OTL 35 Commanding Officer, Headquarters U. S. Army Transportation Research & Engineering Command Fort Eustis, Virginia ATTN: Chief, TechnicalServices Division 36 Commanding General, Ordnance Tank-Automotive Command Detroit Arsenal 28251 Van Dyke Avenue Center Line, Michigan ATTN: Chief, ORDMC-RRS 37 Commanding General, Army Medical Research & Development Command Main Navy Building, Washington 25, D. C. ATTN: Neurophychiatry & Psychophysiology Research Branch 38 Director, U. S. Army Engineer Research & Development Laboratories Fort Belvoir, Virginia ATTN: Chief, Topographic Engineer Department Copy No. Addressee 39-41 Director, U. S. Army Engineer Research & Development Laboratories Fort Belvoir, Virginia ATTN: Chief, Electrical Engineering Department 42 Director, U. S. Army Engineer Research & Development Laboratories Fort Belvoir, Virginia ATTN: Technical Documents Center 43 Commandant, U. S. Army War College Carlisle Barracks, Pennsylvania ATTN: Library 44 Commandant, U. S. Army Command & General Staff College Fort Leavenworth, Kansas ATTN: Archives 45-46 Assistant Commandant, U. S. Army Artillery & Missile School Fort Sill, Oklahoma 47 Assistant Commandant, U. S. Army Air Defense School Fort Belvoir, Virginia ATTN: Combat Development Group 48 Commandant, U. S. Army Aviation School Fort Rucker, Alabama 49 Commanding Officer, U. S. Army Signal Electronic Research Unit P. 0. Box 205 Mountain View, California 50 Office of Naval Operations Department of the Navy, Washington 25, D. C. ATTN: OP-07T 51-53 Office of Naval Research Department of the Navy 17th & Constitution Ave., N. W. Washington 25, D. C. ATTN: Code 463 54 Chief, Bureau of Ships Department of the Navy, Washington 25, D. C. ATTN: Code 312 55-56 Director, U. S. NavalResearch Laboratory Washington 25, D. C. ATTN: Code 2027 57 Commanding Officer, U. S. Navy Ordnance Laboratory Corona, California ATTN: Library 58 Commanding Officer & Director U. S. Navy Electronics Laboratory San Diego 52, California ATTN: Library 49

The University of Michigan Willow Run Laboratories DISTRIBUTION LIST 5 1 June 1960-Effective Date Copy No. Addressee 59 Department of the Air Force Headquarters, USAF Washington 25, D. C. ATTN: Directorate of Requirements 60 Commander, Air Technical Intelligence Center Wright-Patterson Air Force Base, Ohio 61-70 ASTIA (TIPCR) Arlington Hall Station Arlington 12, Virginia 71-75 Commander, Wright Air Development Center Wright-Patterson Air Force Base, Ohio ATTN: WCLROR 76 Commander, Wright Air Development Center Wright-Patterson Air Force Base, Ohio ATTN: WCLDRFV 77 Commander, Wright Air Development Center Wright-Patterson Air Force Base, Ohio ATTN: WCOSI-Library 78 Commander, Rome Air Development Center Griffiss Air Force Base, New York ATTN: RCVSL-1 79 Commander, Rome Air Development Center Griffiss Air Force Base, New York ATTN: RCVH 80-81 Commander, Air Force Cambridge Research Center Laurence G. Hanscom Field Bedford, Massachusetts ATTN: CRES, Stop 36 82-85 Central Intelligence Agency 2430 E. Street, N. W., Washington 25, D. C. ATTN: OCR Mail Room 86-90 National Aeronautics & Space Administration 1520 H. Street, N. W. Washington 25, D. C. 91 U. S. Army Air Defense Human Research Unit Fort Bliss, Texas ATTN: Library 92-93 Combat Surveillance Project Cornell Aeronautical Laboratory, Incorporated Box 168, Arlington 10, Virginia ATTN: Technical Library Copy No. Addressee 94 The RAND Corporation 1700 Main Street Santa Monica, California ATTN: Library 95 Chief, U.S. Army Armor Human Research Unit Fort Knox, Kentucky ATTN: Administrative Assistant 96 Director of Research, U. S. Army Infantry Human Research Unit P. 0. Box 2086, FortBenning, Georgia 97 Chief, U. S. Army Leadership Human Research Unit P. 0. Box 787 Presidio of Monterey, California ATTN: Librarian 98 Chief Scientist, Research &Development Division, Office of the Chief Signal Officer Department of the Army, Washington 25, D. C. 99 Stanford Research Institute Document Center Menlo Park, California ATTN: Acquisitions 100 Operations Research Office The Johns Hopkins University 6935 Arlington Road Bethesda, Maryland, Washington 14, D. C. ATTN: Chief, Intelligence Division 101-102 Cornell Aeronautical Laboratory, Incorporated 4455 Genesee Street, Buffalo 21,NewYork ATTN: Librarian VIA: Bureau of Aeronautics Representative 4455 Genesee Street Buffalo 21, New York 103-104 Control Systems Laboratory University of Illinois Urbana, Illinois ATTN: Librarian VIA: ONR Resident Representative 1209 W. Illinois Street Urbana, Illinois 105-106 Director,Human Resources Research Office The George Washington University P. 0. Box 3596, Washington25, D. C. ATTN: Library 107 Massachusetts Institute of Technology, Research Laboratory of Electronics Cambridge 39, Massachusetts ATTN: Document Room 26-327 50

The University of Michigan Willow Run Laboratories DISTRIBUTION LIST 5 Copy No. Addressee 108 The U. S. Army Aviation HRU P. 0. Box428, Fort Rucker, Alabama 109-110 Visibility Laboratory, Scripps Institution of Oceanography University of California San Diego 52, California 111-113 Bureau of Aeronautics Department of the Navy, Washington 25, D. C. ATTN: RAAV-43 114 Office of Naval Research Department of the Navy 17th & Constitution Ave., N. W. Washington 25, D. C. ATTN: Code 461 115 Director, Electronic Defense Group U of M Research Institute The University of Michigan Ann Arbor, Michigan 1 June 1960- Effective Date Copy No. Addressee 116-118 Assistant Commandant U. S. Army Air Defense School Fort Bliss, Texas 119-120 Mitre Corporation P. 0. Box 208 Lexington 73, Massachusetts VIA: Commander Air Defense Systems Integration Division United States Air Force Laurence G. Hanscom Field Lexington 73, Massachusetts 121 U. S. Continental Army Command Liaison Officer Project MICHIGAN, Willow Run Laboratories Ypsilanti, Michigan 122 Commanding Officer U. S. Army Liaison Group Project MICHIGAN, Willow Run Laboratories Ypsilanti, Michigan 51

I AD Div. 15/1 Willow Run Laboratories, U. of Michigan, Ann Arbor A STOCHASTIC METHOD OF SOLVING PARTIAL DIFFERENTIAL EQUATIONS USING AN ELECTRONIC ANALOG COMPUTER (Unclassified title) by Kuei Chuang, Louis F. Kazda, and Thomas Windeknecht. Report of Project MICHIGAN. Jun 60. 48 p. incl. illus., 6 refs. (Report no. 2900-91-T) (Contract DA-36-039 SC-78801) Unclassified report This report presents an analog-computer method for solving a general class of two-dimensional boundary-value problems in which only a small number of operational amplifiers is required. The method simulates the boundary-value problem under consideration as a stochastic model by a technique commonly designated the "Monte Carlo method. " In this case, the stochastic method becomes a tool for the numerical analysis of the partial differential equation. Although the Monte Carlo method of digital-computer solution of partial differential equations has been extensively studied, its application to the analog computer is believed to be new. The boundary-value problems for which the stochastic method is applicable belong to a family of generalized Dirichlet problems. (over) + UNCLASSIFIED 1. Mathematical computers — Analog computers 2. Mathematical computersApplications I. Title: Project MICHIGAN II. Chuang, Kuei; Kadza, Louis F.; Windeknecht, Thomas III. U. S. Army Signal Corps IV. Contract DA-36-039 SC-78801 Armed Services Technical Information Agency UNCLASSIFIED UNCLASSIFIED 1. Mathematical computersAnalog computers 2. Mathematical computers — Applications I. Title: Project MICHIGAN II. Chuang, Kuei; Kadza, Louis F.; Windeknecht, Thomas III. U. S. Army Signal Corps IV. Contract DA-36-039 SC-78801 Armed Services Technical Information Agency UNCLASSIFIED AD Div. 15/1 Willow Run Laboratories, U. of Michigan, Ann Arbor A STOCHASTIC METHOD OF SOLVING PARTIAL DIFFERENTIAL EQUATIONS USING AN ELECTRONIC ANALOG COMPUTER (Unclassified title) by Kuei Chuang, Louis F. Kazda, and Thomas Windeknecht. Report of Project MICHIGAN. Jun 60. 48 p. incl. illus., 6 refs. (Report no. 2900-91-T) (Contract DA-36-039 SC-78801) Unclassified report This report presents an analog-computer method for solving a general class of two-dimensional boundary-value problems in which only a small number of operational amplifiers is required. The method simulates the boundary-value problem under consideration as a stochastic model by a technique commonly designated the "Monte Carlo method. " In this case, the stochastic method becomes a tool for the numerical analysis of the partial differential equation. Although the Monte Carlo method of digital-computer solution of partial differential equations has been extensively studied, its application to the analog computer is believed to be new. The boundary-value problems for which the stochastic method is applicable belong to a family of generalized Dirichlet problems. (over) AD Div. 15/1 Willow Run Laboratories, U. of Michigan, Ann Arbor A STOCHASTIC METHOD OF SOLVING PARTIAL DIFFERENTIAL EQUATIONS USING AN ELECTRONIC ANALOG COMPUTER (Unclassified title) by Kuei Chuang, Louis F. Kazda, and Thomas Windeknecht. Report of Project MICHIGAN. Jun 60. 48 p. incl. illus., 6 refs. (Report no. 2900-91-T) (Contract DA-36-039 SC-78801) Unclassified report This report presents an analog-computer method for solving a general class of two-dimensional boundary-value problems in which only a small number of operational amplifiers is required. The method simulates the boundary-value problem under consideration as a stochastic model by a technique commonly designated the "Monte Carlo method. " In this case, the stochastic method becomes a tool for the numerical analysis of the partial differential equation. Although the Monte Carlo method of digital-computer solution of partial differential equations has been extensively studied, its application to the analog computer is believed to be new. The boundary-value problems for which the stochastic method is applicable belong to a family of generalized Dirichlet problems. (over) UNCLASSIFIED 1. Mathematical computers — Analog computers 2. Mathematical computers — Applications I. Title: Project MICHIGAN II. Chuang, Kuei; Kadza, Louis F.; Windeknecht, Thomas III. U. S. Army Signal Corps IV. Contract DA-36-039 SC-78801 Armed Services Technical Information Agency UNCLASSIFIED UNCLASSIFIED 1. Mathematical computersAnalog computers 2. Mathematical computers — Applications I. Title: Project MICHIGAN II. Chuang, Kuei; Kadza, Louis F.; Windeknecht, Thomas III. U. S. Army Signal Corps IV. Contract DA-36-039 SC-78801 Armed Services Technical Information Agency UNCLASSIFIED AD Div. 15/1 Willow Run Laboratories, U. of Michigan, Ann Arbor A STOCHASTIC METHOD OF SOLVING PARTIAL DIFFERENTIAL EQUATIONS USING AN ELECTRONIC ANALOG COMPUTER (Unclassified title) by Kuei Chuang, Louis F. Kazda, and Thomas Windeknecht. Report of Project MICHIGAN. Jun 60. 48 p. incl. illus., 6 refs. (Report no. 2900-91-T) (Contract DA-36-039 SC-78801) Unclassified report This report presents an analog-computer method for solving a general class of two-dimensional boundary-value problems in which only a small number of operational amplifiers is required. The method simulates the boundary-value problem under consideration as a stochastic model by a technique commonly designated the "Monte Carlo method. " In this case, the stochastic method becomes a tool for the numerical analysis of the partial differential equation. Although the Monte Carlo method of digital-computer solution of partial differential equations has been extensively studied, its application to the analog computer is believed to be new. The boundary-value problems for which the stochastic method is applicable belong to a family of generalized Dirichlet problems. (over)

AD UNCLASSIFIED The analog computer, using the stochastic method, solves the UNITERMS boundary-value problem at a desired point. This solution, which Partal differential equation can be shown to be a statistically convergent quantity, appears as Monte Carlo method the mean value of the summation of discrete boundary values that Analog computers Analog computers are randomly generated by the computing system. Dirichlet A number of one- and two-dimensional boundary-value problems Stochastic of the type described, whose analytical solutions are known, have Boundary-value problem been experimentally solved. A comparison of the experimental and analytical solutions corroborates the validity of this method. AD UNCLASSIFIED The analog computer, using the stochastic method, solves the UNITERMS boundary-value problem at a desired point. This solution, which Paal dfferential equation can be shown to be a statistically convergent quantity, appears as M e Carlo m d the mean value of the summation of discrete boundary values that Anao com e are randomly generated by the computing system. Dirichlet A number of one- and two-dimensional boundary-value problems Stochastic of the type described, whose analytical solutions are known, have Boundary-value problem been experimentally solved. A comparison of the experimental and analytical solutions corroborates the validity of this method. AD The analog computer, using the stochastic method, solves the boundary-value problem at a desired point. This solution, which can be shown to be a statistically convergent quantity, appears as the mean value of the summation of discrete boundary values that are randomly generated by the computing system. A number of one- and two-dimensional boundary-value problems of the type described, whose analytical solutions are known, have been experimentally solved. A comparison of the experimental and analytical solutions corroborates the validity of this method. UNCLASSIFIED UNITERMS Partial differential equation Monte Carlo method Analog computers Dirichlet Stochastic Boundary-value problem AD The analog computer, using the stochastic method, solves the boundary-value problem at a desired point. This solution, which can be shown to be a statistically convergent quantity, appears as the mean value of the summation of discrete boundary values that are randomly generated by the computing system. A number of one- and two-dimensional boundary-value problems of the type described, whose analytical solutions are known, have been experimentally solved. A comparison of the experimental and analytical solutions corroborates the validity of this method. UNCLASSIFIED UNITERMS Partial differential equation Monte Carlo method Analog computers Dirichlet Stochastic Boundary-value problem

I AD Div. 15/1 Willow Run Laboratories, U. of Michigan, Ann Arbor A STOCHASTIC METHOD OF SOLVING PARTIAL DIFFERENTIAL EQUATIONS USING AN ELECTRONIC ANALOG COMPUTER (Unclassified title) by Kuei Chuang, Louis F. Kazda, and Thomas Windeknecht. Report of Project MICHIGAN. Jun 60. 48 p. incl. illus., 6 refs. (Report no. 2900-91-T) (Contract DA-36-039 SC-78801) Unclassified report This report presents an analog-computer method for solving a general class of two-dimensional boundary-value problems in which only a small number of operational amplifiers is required. The method simulates the boundary-value problem under consideration as a stochastic model by a technique commonly designated the "Monte Carlo method. " In this case, the stochastic method becomes a tool for the numerical analysis of the partial differential equation. Although the Monte Carlo method of digital-computer solution of partial differential equations has been extensively studied, its application to the analog computer is believed to be new. The boundary-value problems for which the stochastic method is applicable belong to a family of generalized Dirichlet problems. (over) AD Div. 15/1 Willow Run Laboratories, U. of Michigan, Ann Arbor A STOCHASTIC METHOD OF SOLVING PARTIAL DIFFERENTIAL EQUATIONS USING AN ELECTRONIC ANALOG COMPUTER (Unclassified title) by Kuei Chuang, Louis F. Kazda, and Thomas Windeknecht. Report of Project MICHIGAN. Jun 60. 48 p. incl. illus., 6 refs. (Report no. 2900-91-T) (Contract DA-36-039 SC-78801) Unclassified report This report presents an analog-computer method for solving a general class of two-dimensional boundary-value problems in which only a small number of operational amplifiers is required. The method simulates the boundary-value problem under consideration as a stochastic model by a technique commonly designated the "Monte Carlo method. " In this case, the stochastic method becomes a tool for the numerical analysis of the partial differential equation. Although the Monte Carlo method of digital-computer solution of partial differential equations has been extensively studied, its application to the analog computer is believed to be new. The boundary-value problems for which the stochastic method is applicable belong to a family of generalized Dirichlet problems. (over) UNCLASSIFIED 1. Mathematical computers — Analog computers 2. Mathematical computers — Applications I. Title: Project MICHIGAN H. Chuang, Kuei; Kadza, Louis F.; Windeknecht, Thomas III. U. S. Army Signal Corps IV. Contract DA-36-039 SC-78801 Armed Services Technical Information Agency UNCLASSIFIED UNCLASSIFIED 1. Mathematical computersAnalog computers 2. Mathematical computers — Applications I. Title: Project MICHIGAN II. Chuang, Kuei; Kadza, Louis F.; Windeknecht, Thomas III. U. S. Army Signal Corps IV. Contract DA-36-039 SC-78801 Armed Services Technical Information Agency UNCLASSIFIED AD Div. 15/1 Willow Run Laboratories, U. of Michigan, Ann Arbor A STOCHASTIC METHOD OF SOLVING PARTIAL DIFFERENTIAL EQUATIONS USING AN ELECTRONIC ANALOG COMPUTER (Unclassified title) by Kuei Chuang, Louis F. Kazda, and Thomas Windeknecht. Report of Project MICHIGAN. Jun 60. 48 p. incl. illus., 6 refs. (Report no. 2900-91-T) (Contract DA-36-039 SC-78801) Unclassified report This report presents an analog-computer method for solving a general class of two-dimensional boundary-value problems in which only a small number of operational amplifiers is required. The method simulates the boundary-value problem under consideration as a stochastic model by a technique commonly designated the "Monte Carlo method. " In this case, the stochastic method becomes a tool for the numerical analysis of the partial differential equation. Although the Monte Carlo method of digital-computer solution of partial differential equations has been extensively studied, its application to the analog computer is believed to be new. The boundary-value problems for which the stochastic method is applicable belong to a family of generalized Dirichlet problems. (over) AD Div. 15/1 Willow Run Laboratories, U. of Michigan, Ann Arbor A STOCHASTIC METHOD OF SOLVING PARTIAL DIFFERENTIAL EQUATIONS USING AN ELECTRONIC ANALOG COMPUTER (Unclassified title) by Kuei Chuang, Louis F. Kazda, and Thomas Windeknecht. Report of Project MICHIGAN. Jun 60. 48 p. incl. illus., 6 refs. (Report no. 2900-91-T) (Contract DA-36-039 SC-78801) Unclassified report This report presents an analog-computer method for solving a general class of two-dimensional boundary-value problems in which only a small number of operational amplifiers is required. The method simulates the boundary-value problem under consideration as a stochastic model by a technique commonly designated the "Monte Carlo method. " In this case, the stochastic method becomes a tool for the numerical analysis of the partial differential equation. Although the Monte Carlo method of digital-computer solution of partial differential equations has been extensively studied, its application to the analog computer is believed to be new. The boundary-value problems for which the stochastic method is applicable belong to a family of generalized Dirichlet problems. (over) UNCLASSIFIED 1. Mathematical computers — Analog computers 2. Mathematical computers — Applications I. Title: Project MICHIGAN II. Chuang, Kuei; Kadza, Louis F.; Windeknecht, Thomas III. U. S. Army Signal Corps IV. Contract DA-36-039 SC-78801 Armed Services Technical Information Agency UNCLASSIFIED UNCLASSIFIED 1. Mathematical computersAnalog computers 2. Mathematical computersApplications I. Title: Project MICHIGAN II. Chuang, Kuei; Kadza, Louis F.; Windeknecht, Thomas III. U. S. Army Signal Corps IV. Contract DA-36-039 SC-78801 Armed Services Technical Information Agency UNCLASSIFIED

AD UNCLASSIFIED The analog computer, using the stochastic method, solves the UNITERMS boundary-value problem at a desired point. This solution, which Partial differential equation can be shown to be a statistically convergent quantity, appears as Monte Carlo method the mean value of the summation of discrete boundary values that Analog computers are randomly generated by the computing system. Dirichlet A number of one- and two-dimensional boundary-value problems Stochastic of the type described, whose analytical solutions are known, have Boundary-value problem been experimentally solved. A comparison of the experimental and analytical solutions corroborates the validity of this method. AD UNCLASSIFIED The analog computer, using the stochastic method, solves the UNITERMS boundary-value problem at a desired point. This solution, which Partial differential equation can be shown to be a statistically convergent quantity, appears as Monte Carlo method the mean value of the summation of discrete boundary values that Analog computers are randomly generated by the computing system. Dirichlet A number of one- and two-dimensional boundary-value problems Stochastic of the type described, whose analytical solutions are known, have Boundary-value problem been experimentally solved. A comparison of the experimental and analytical solutions corroborates the validity of this method. AD UNCLASSIFIED The analog computer, using the stochastic method, solves the UNITERMS boundary-value problem at a desired point. This solution, which P l d Partial differential equation can be shown to be a statistically convergent quantity, appears as Monte Crlo ethod the mean value of the summation of discrete boundary values that Analo omter are randomly generated by the computing system. Dirichlet A number of one- and two-dimensional boundary-value problems Stochastic of the type described, whose analytical solutions are known, have Boundary-value problem been experimentally solved. A comparison of the experimental and analytical solutions corroborates the validity of this method. AD UNCLASSIFIED The analog computer, using the stochastic method, solves the UNITERMS boundary-value problem at a desired point. This solution, which Partial differential equation can be shown to be a statistically convergent quantity, appears as Monte Carlo method the mean value of the summation of discrete boundary values that Analog computers are randomly generated by the computing system. Dirichlet A number of one- and two-dimensional boundary-value problems Stochastic of the type described, whose analytical solutions are known, have Boundary-value problem been experimentally solved. A comparison of the experimental and analytical solutions corroborates the validity of this method. Co - c CO -z 01 -= o - X 0 1 -G) Ul;-