PLOD: Interative solver for ordinary differential equations

302

PLOD: I n t e r a c t i v e S o l v e r for O r d i n a r y Differential E q u a t i o n s David K. Kahaner and Selden L. Stewart National Institute of Standards and Technology Center for Computing and Applied Mathematics Technology Building, Room A151 Gaithe~sburg, MD 20899 :-.,~vid Barnett University of the District of Columbia Computer Science Department Washington, DC A b s t r a c t : PLOD is an acronym for PLotted solutions of Differential Equations. PLOD can solve up to 25 ordinary differential equations with 10 embedded parameters. It is entirely interactive, requiring little progrAmmlcg effort. PLOD uses current numerical methods, as implemented in the DDRIV package, to perform the integration. This includes the ability to solve "stiff" equations. The current version of PLOD has been implemented on the IBM Personal Computer family, such as the XT, AT and PS/2. PLOD is in the public domain.

L Introduction In this paper we consider the practical issues involved in studying a physical model which is described by a set of ordinary differential equations (ODEs) with prescribed initial conditions; an initial value problem. Many such problems involve a doze.n or so differential equations and a similar number of parameters. These problems occur frequently in studying dynamics or in the physical modelling of such phenomena as Josephson junctions [I]. They can be among the most taxing to solve numerically because of potential stiffness and instability with respect to initial conditions. Few models can be solved explicitly in closed form; most cannot. Some are easy for almost any numerical method. These often occur in university situations as examples for students. There are several interactive ODE solvers which address this class of problem, for instance, Phaser [2], but because they use relatively primitive numerics, such solvers are unsuitable for more difficult cases. Simulation packages may incorporate good ODE solvers, but in a compendium published in uSimuiation" [3], none of those listed were in the public domain. A few of the higher priced packages make use of well designed integrators. In some packages there is a large choice of basic integrator, and virtually all packages allow for some choice. A scientist at a university or research laboratory who purchases a product in the $400-$700 range for occasional use is likely to get one that does not incorporate state of the art integration techniques. Common integrators are classical explicit Runge-Kutta or RungeKutta Fehlberg, which are ineffective for the solution of stiff differential equations. On the other hand, the more powerful (and costly) simulation packages are usually large, and can be difficult to use. These may have extensive capabilities that are crucial in some applications, hut often go beyond the range of needs of many users. Most have detailed reference manuals, often in several volumes, and can have their own language syntax to be mastered. A few are batch oriented. By that we mean that the user describes the problem by means of a file which is then input to the simulator. It may not be possible to make simple changes in some parameter (unless this is known in advance) without rerunning the job. Many packages have the capability of graphical output but some still do not. 0096-3003189/$00.00

303

There seemed to us to be a need for software which addressed the more casual use., but which was of very high scientific quality. Furthermore, scientists who actually solve differential equations are likely to have a powerful computer on their desks, usually a 32 bit Intel 80386/80387 based "PC." The PLOD package (PLotted solutions of Ordinary Differential equations) has been designed with the following criteria: (1) The physical problem should be able to be described in terms of 25 (or fewer) ordinary differential equations with up to 10 parameters. (2) The problem should be such that integration can be done in real time. PLOD is entirely interactive. Item number (1) addresses a large fraction of the problems which arise in "applications '~, and many prototype problems in dynamics. However, PLOD is inappropriate for studying, say, the complete suspension system of an automobile. Item (2) can also exclude small problems, such as those which might occur in "chaotic" [4] systems where there are many abrupt changes in the dependent variables, or where the user wants to integrate for long periods of simulated time. It also excludes problems in which the user wants to repeat a computation many times (automatically) to compute an average as some parameter is methodically changed. Thus PLOD is designed to deal with "what if*' type problems, experimental studies, and perhaps as a prototype for more ambitious simulations [4]. (3) PLOD should be easy enough to use that any scientist can get it to run without a manual. A text file containing installation and tutorial information is included in the PLOD distribution. It is useful to re.~d, but once the package has been installed, even first time users can function without it. Inadvertent keyboard blunders should not prevent the package from working. It is never completely possible to anticipate all combinations of sltcations, but error recovery must be a major design goal. Thus PLOD is well suited for use in a university situation where students will explore the effects of changing a model. (4) The nume~cal methods and their implemevtation should be current. The underlying integrator, DDRIV [5], has been in use over a number of years at numerous sites It uses modem (but by now, quite conventional) methods including careful implementation of variable step and ~arlable order Adams and Gear formulas [0], which give it the ability to solve stiff equations. Other integrators would be equally suitable, and could he exchanged for DDRIV, but we think of PLOD as a tool for solving problems rather than an environment for experimenting with different solvers. (5) Rapid, flexible and attractive graphics should be an integral part of the package. On the other hand, publication qughty graphics is not included as a design criteria. (6) As much as possible, PLOD should be portable. Th.~t is, special charar.teristics of a particular computer should not be used l~r2ess these are unavoidable. Beca,me PLOD utilizes graphics, it is not entirely portable, but non portable section~ should be kept to a minimum and isolated as much as possible. (7) All aspects of PLOD should be in the public domain. Thus the package is completely free to use, copy and distribute. This paper describes the major characteristics of PLOD from the perspective of a potential user.

304

H. Overview of PLOD PLOD solves =

y

(a) =

(1) = Y--

The solution is aesired on an interval [a, b]. In most problems of physical interest the equations also involve constants called physical parameters whose values determine the character of the solution. Equation (1) does not show these explicitly. PLOD is designed to be used in two steps, called PLOD0 and PLOD1. We think of these steps as ~modd building ~ and "model analysis" respectively. Both are interactive. A small "batch" program, called PLOD, which the user executes by typing its aame, runs PLOD0 and optionally will do the processing to pass to PLOD1. Rlmnlng PLOD0, the user is queried for a description of the model, i.e., the equations. No integration occurs during this step and1 in fact, numerical values are not requested for any of the variables or parameters during the model building phase. PLOD0 is a preprocessor; it~ cutp~.lt is a Fortran program having the user's equations in a form which will enable them to be integrated numerically. It also has a mechanism to communicate the names and other information to PLOD1. In order to continue, the user takes the Fortran program which was output from PLOD0, compiles it and then links to several already precompiled modules which are supplied with the package. The process of compiling and linking can be done automatically vi~ the batch program, PLOD. Alternatively, PLOD0 also gives instructions for manual compiling and linking. The result is another executable program which the user can then execute. It is the running of this latter program which we call PLOD1. During PLOD] the user is asked for values of the parameters, initial values for the dependent variables, and the interval on which the integration is to take place. After the integration the user can request various information about any combinations of variables and parameters, either through plots, Poincare plo~s, or listings. It is also possible to change the parameter values, initial conditions, integration interval, etc., and experiment with the problem, examining various results under varying conditions. A more sophisticated user can also in:eractively alter many of DDRIV's parameters, including the integration and iteration methods, error tolerances, and so forth. During PLOD] the user cannot alter the functional form of the mode] by adding terms, or adding or removing equations. To make those changes requires returning to PLOD0. HI. W h a t D o Y o u N e e d t o R u n P L O D ? The first version of PLOD was mainframe oriented, and tested on a Sperry ]]00/80 at the National Bureau of Standards (NBS) as well as a Digital Equipment Vax 11/780 at Los Alamos National Laboratory. It was written entirely in Fortran 77 except that it used DISSPLA [7] to generate all 9f its graphics. The current version has been developed for "PC" type computers with at least 3fi0K of memory. A typical executable file generated by PLOD1 averages about 320K for 5 to 10 equations. In practice PLOD works best on high-end computers. For example, to integrate through one full cycle of a model of a chemical oscillator (see Figure ]) [8, p.50] requires about 4 seconds on an IBM PS/2 model 80, 20 MHz. Graphics is an integral par~ of PLOD, which supports the IBM Color Graphics Adaptor (CGA), IBM Enhanced Graphics Adaptor and Enhanced Color Display (EGA),

305

IBM Virtual Device Array (VGA) and the Hercules Graphics Card [9] (HGA). PLOD does not burden the user with the need to specify which graphics card was in use, but will run automatically using the hi~hest resolution, on an CGA, E t A , VGA or HGA. Because the IIGA does not have any color capabilities and because some IBM compatible computers do not support color, PLOD does not utilize color in any way. PLOD is also likely to work on any computer that is true]y "PC ~' compatible. To obtain hard copy output of the plots, a graphics printer is required. PLOD supports most standaxd IBM or Epson model~. PLOD can produce "screen quality" plots on the printer. To obtain publication quality plots a user would normally save the integration data and postprocess it with locally available software. The figures shown in this paper were generated in this manner. The latebt version of PLOD requires the u ~ of Lahey Computer Systems (F77L) Fortran compiler [10], version 2.2 and a compatible version of the IBM Linker, version 2.3 or above. A version of PLOD is also available to use with Ryan Mc-Farland Fortran [11] (RMFort] version 2.4 and Phoenix Software's PLINK86 linker [12]. We chose lo implement in Fortran rather than C or Pascal because of familiarity (by us and our user community) and because the best ODE solvers are already in the public domain, in Fortran. A recent book [8] lists more thzm 50 pages of references to research papers and software. None of the software references mention any language except Fortran° PLOD0 (the preprocessor) is written in Pascal; most of PLOD1 (the integration and analysis step) is written in Fortran 77. In fac~. the PLOD1 source has been a good test for Fortran vendors. The sections dealing directly with graphics have been written ~n a combination of C and Intel 8088 Assembler. A version of PLOD is also available which does not require either a compiler or a linker. This avoids some installation difficulties and might be suitable for a classroom environment. On the other hand it r,:qvJre~ that the differential equations be interpreted, resulting in much slower integration. For research applications we feel that the compiled version is more appropriate. IV. P L O D O E x a m p l e : P i l o t Seat E j e c t i o n S t u d y To illustrate the basic approach we present below a sample sess:on in the use of PLOD, taking as an example a simplified model of pilot ejection from an aircraSt. We assume that the pilot is strapped into a seat in an aircraft which is flying horizontally at constant velocity Va. At ~ = 0 the seat-pilot combination is forced up and back along rails in the cockpit until the unit leaves the aircraft and subsequently follows a ballistic trajectory. It is important to discover if the initial vertical velocity of the seat :s sufficient for it to clear the vertical stabilizer. If we place an z y coordinate system in the cockr". near the base of the seat then the subsequent motion can be described by the equations

[13], z' = v cos e - V~

Y'

= vsinO

D = pCdSv2/2

v'={

0 -

Dim - 9 sin 0

0'={ 0 -(9/~)cos0

ifO<_y W if 0 < Y < ya; if y > ya

306

where z and 9 are the position of the seat relative to the aircraft, v is the magnitude of the seat's velocity vector with respect to ground coordinates, and 0 is its angle to the horizontal. There are also a number of other indeterminates, D, Va, etc. Some of these are inserted for our convenience, such as D, or g. The former is the drag on the seat-vilot combination and is entirely defined in terms of other variables. On the other hand 9 which is thv .~.~:eleration of gravity is a fixed constant. Others, such as Va,p, etc., are pawmeters in the sense that we may want to alter their values to see the effect on the motion. These parameters are: Va: the horizontal speed of the aircraft, assumed constant, yQ: the position above the horizontal axis at which the seat is assumed fully ejected, Cd: the drag coefficient of the seat, S: the average cross-sectional area of the seat presented to the air flow, p: the air density, and m: the mass of the pilot-seat combination. Thus there are n = 4 ODEs. There are in adn'don, 6 parameters Ya,p, Cd, S,y~ and m. We do not think of g as a parameter since its value is never expected to change. Similarly, D is also not considered a parameter because its value can al 7ays be determined explicitly by the values of the other variables and parameters. Reasonable initial values are • (0) = 0

v(0) = o ,,(o) = v/(V

- v, smo,)2 + (v, cos

(2)

[v.coso. ] 0(0) = tan -1 [V~ - Ve sin0e] ' where Ve is the speed of ejection and 0~ is the angle of eject_;on measured from the vertical (the rails on which the seat travels are tilted slightly backwards). PLOD comes with two tutorial fdes which provide the details, but v'e believe that, except for reading the installation instructions it is possible to proceed without them. Installation involves copying all the fdes onto a directory and altering one llne of the file which is used in the automatic linking step. Subsequent to that, a user can b e ~ n by typing PLOD. In what follows we denote user responses in boldface. Our comments, which obviously do not appear in the program, are included in brackets [ ] where appropriate. PLOD PLOD here for solving differential equations... Help? (Y/N) iN] n [A user responding affirmatively to this question is given some brief information about the package, including a llst of the necessary hardware and software, and referred to the tutorials. All questions have default responses whi:h are indicated in the brackets, as {N] above, and is obtained by entering a carriage return.] Enter FILE name on which P L O D will store the Fortran program it generates. PILOT lit is possible to manipulate a model file that was previously created, or to create a model interactively. In either case PLOD drops the user into a "model editor." First time users proceed interacti~ely L.,L ,,bu,dly bare ,heir model for later alteration. Figure 2 shows the model editor after the user has entered the equations. These mostly follow ordinary Fortran syntax, except that the primed quantities, such as z~, etc., may appear either on the left of assignment statements, or subsequently, anywhere else.]

307

To MOVE CURSOR u s e A r r o v s , PSUp, PgDn, Home, End To EXIT k CONTINUE u s e E s c , To F.EASE CHARACTER ABOVE CURSOR u s e DEL

To DELETE FROM CURSOR TO LIXE END i s e C t r ; and T ( t o g e t h e r ) Note: Independent v a r i a b l e i s T [ ........

1 .........

2 .........

3,

.4,

5,

6

co~ments begin with p e r c e n t s i g n ~. and i g n o r e r e s t o f l i n e . Pilot seat ejection study Y. x , y : Coordinate system in c o c k p i t

Y. N o t e :

v :

Magnitude of seat's velocity vector vrt

ground ~. t h e t a : A n g l e of s e a t t o h o r i z o n t a l . ~, Va : Uorizontal aircraft speed ~, ya : Position above horizontal at vhich seat i s fully e j e c t e d ~, Cd : Drag coefficient of seat ~, S : Average cross section of seat wrt air flow

Y. rho : m :

Air d e n s i t y Mass o f p t l o t + s e a ~

x ~=v*cos(theta)-va y ' =v* s i n ( t h e t a) Y, gravity if (y. it.ya) then v2~O

g=32.2

t h e t a ~= 0 else D=. 5*rho*Cd*S*v*v v i=-D/m-g*sin(thet a)

theta' =-g,cos (theta) / endif

F i g u r e 2: P L O D ' s M o d e l E d i t o r [After exiting the editor and checking model syntax PLOD will cause the program which has been generated to he compiled, linked to the necessary libraries, and finally to begin execution of PILOT.EXE, which embodies the PLOD1 step.] V. P L O D 1 E x a m p l e : P i l o t S e a t E j e c t i o n S t u d y PILOT is menu driven, i.e., the user is asked to select from one of a number of possible choices displayed on the screen. So every interactive session will differ in the order in which operations are performed, hut typically, users will input their data, integrate, plot and or llst the integration results, modify some aspect of the problem, reexamine the output, etc. The current user has one important parameter change to investigate, via, Va = 900, and V= = 500 (feet/second). The key questiov for the designer is whether the ejected seat c]ears the rear stabilizer for this range of V=. Reasonable values for the other parameters are y= = 4 ft, m = 7 slugs, p : 2,3769.10 -s slugs/ft s, C¢~ = 1, S = 10 ft 2, g = 32.2 ft/sec 2, along with initial values given in Equation (2), where Ve = 40 ft/sec is the ejection velocity of the seat, and 8e = 15/57.3 (radians) is the ejection angle.

308

[PLOD1 begins by displaying the ODEs (omitted here for brevity), notes some simple defaults, and oilers users an opportunity to ret.d values that might have been saved from an earlier ses~inn.] Defaults: Integration interval= [0,1], In:~tia| value(s)=l.0, Parameter(s)=0.0 Want to read problem input values from a file previously written? (Y/N) IN] n Alter any (internal) integration/plotting defaults? (Y/N) IN] n [_Allfirst time asers will respond N to this query, but entering Y gives an opportunity to modify items from an aExpert" menu, which is shown below. PLOD Expert Menu Title: Pilot Seat Ejection Study E;ps: accuracy of integration: 1.000000F_~04 Max integration steps before warning: 400 Plotting tolerance: 1.000000E-03 percent Frequency of saving points: 1 Connect points in plot: yes Iteration method: Newton's method N:Integration method: Gear's method O:Maxhnum integration order: 5 S:Maximum integration step size: Defaulted Z:Problem Zero: 1.110000E-15 B:Variables kept within bounds: Defaulted L:Cycles for non-positive log plots: 4 Add dots to plots: no Data about integration G:Digits in displeys: 4 Newton's method requires finite difference approximations to Jacohian matrix. Functional iteration does not and is faster, but is unsuitable for stiff problems. Pressing space bar toggles from one to the other. Select w i t h arrows or A i t + f i r s t c h a r a c t e r ; Ese to exit Using this menu, it is posfible to alter the integration accuracy, switch between Gear's stiffly stable integration method and an Adams Moulton method, or select either Newton's method or functional iteration for solving the nonlinear implicit equations. It is also possible to adjust some details such as the frequency of saving integration output, number of cycles for log plots, etc. The N response generates a Modify menu, below, in which default parameter and initial values can be altered. ] OPTIONS: (Enter one letter to) Sl:ow equations, initial values, etc alter initial Values at 0.00 alter stopping conditions on... Independent variable i.e., alter interval [0.000, 1.000] Dependent variable(s) e.g., when X=3.14159 (off) alter Parameters Va, etc. Toggle into quiet mode integrate l:[estart the problem

309

[Since neither the initial conditions nor parameters have yet been set, now is the time for this. For brevity, we omit the details of setting the parameter values and assume t h e user has gone on to set the initial values.] V [PLOD responds with a list of variables with their current values, and asks for the names of those to be altered (responses omitted).] Enter initial value for x or < R E T U R N > 0.0 y or < R E T U R N > 0.0 v or < R E T U R N > sqrt((Va-40*sin(15/57.8))**2+(40*eos(15/57.3))**2) theta or < R E T U R N > atan(40*eos(15/ST.8)/(Va-40*sln(15]5T.3))) [In all cases where PLOD expects numerical input, an expression can also be entered involving any of the variables or parameters. PLOD's parser will catch syntax errors, square roots of negative numbers, etc.] [The modify menu is displayed again (omitted). Normal procedure now would he to enter G for integrate, but instead our user types D to specify that the integration should stop when one of the depende,~f variables reaches a particular value. In the context of the sample problem, we want to terminate the computation when the seat moves past the end of the aircraft. Since the initial condition will set z -- 0, we want to terminate when z < -60.] D

Want to specify any stopping criteria based on dependent variables? ( Y / N / ? ) y [A response of ? provides advice.] Enter ? or E x p r l = E x p r 2 : x - - - 6 0 [The user has requested that the integration proceed u ~itil either the right endpo'mt of the integration interval has been reached or • + 60 changes sign.] [The Modify menu appears again (omitted), and the user types Gi to integrate. Details of the specific method used to integrate, accuracy of integration, etc., are unnecessary to specify. (Users who have ]eanlngs in these directions can make changes via the Expert menu.) During the integration, intermediate output is displayed one screenful at a time, to allow verification that things are behaving properly. This output can he shut down by toggling Quiet mode-whlch most users do, but is invaluable at the outset. For nonlinear ODEs, seemingly trivial blunders in problem specification can lead to serious difl~.culties, and this is one way to verify that the solutions are going in the correct direction. Additionally, we have ~ried to design-in crash protection. For example, PLOD will take at most 400 internal steps before halting and asking "what-next?" (400 can be changed in the Expert menu.) Even in Quiet mode, by pressing the space bar, a user can intt :upt the calculation, and if all is well, continue. Finally, if the integration is not making progress (perhaps because of the presence of a singularity) PLOD will halt and make some suggestions. We list below a portton of the first screen-full that our examp]e generates. The number of digits displayed (and consequently the number of variables) is ~..ijust able. ]

310

t

x

y

0.000000E+00 0.596046E-10 0.315910E-06 0.631759E-06 0.947609E-06 0.947877E-06 0.232271E-05 0.369754E-05 0.507236E-05 0.507380E-05 0.169902E-04 0.289066F_,-04 0.408299E-04 0.544094E-02 0.108411E-01 0.162412E-01 0.216413E-01 0.324415E-01 0.432418E-01

0.000000E+09 0.000000E+00 -0.617029E-09 0.230296E-08 -0.327031E-05 0.122059F,-04 -0.654000F-~05 0.244094E-04 -0.980969F,-05 0.366130E-04 -0.981146E-05 0.366225F,-04 -0.240437F_~04 0.897421E-04 -0.382760E-04 0.142862E-03 -0.525083E-04 0.195981E-03 -0.525265E-04 0.196036E-03 -0.175885E-03 0.656452E-03 -0.299244E-03 0.111687F_~02 -0A22603F,-03 0.157728E-02 -0.563249F_~01 0.210223E+00 -0.112227E-]-00 0.418869E+00 -0.168129E-F00 0.627514E+00 -0.22a-032E+00 0.836160E+00 -0.335836E+00 0.125345E+01 -0.447640E-t-00 0.167074E+01

v

theta

0.8905E+03 0.8905E+03 0.8905E+03 0.8905E+03 0.8905E+03 0.8905E+03 0.8905E+03 0.8905E+03 0.8905E+03 0.8905E+03 0.8905E+03 0.8905E+03 0.8905E+03 0.8905E+03 0.8905E+03 0.8905E+03 0.8905E+03 0.8905E+03 0.8905E+03

0.434E-01 0.434E-01 0.434E-01 0.434F,-01 0.434E-01 0.434E-01 0.434E-01 0.434E-01 0.434E-01 0.434E-01 0.434E-01 0.434E-01 0.434E-01 0.434E-01 0.434E-01 0.434F,-01 0.434E-01 0.434E-01 0.434E-01

Hit to continue integration, S to stop, or Q to continue quietly. Q Integrating ... (press space bar to halt) [The vaxiabiy spaced time points displayed are those which are output from the integrator. ] J

[Our user has gone into Quiet mode, so his te.-~:nlnal i~ silent until] Integration complete to t=0.433686 because x---60.00000 At this point x = -60.0000000

y= 12.80253124 v= 593.60034180 theta--2.84800231E-02 This stopping criteria remains active. Integrate further toward t=l.0, or Stop integration, or Toggle quiet mode and integrate further, or store Point and continue quietly? (I/S/T/P) [S] s [The 'P' option permits the generation of Poincare plots [14]. The "Main" menu now appears.] Options: (Enter one letter) Show equations, inhial values... Graphics Full

Brief (no legends...) List results on screen Modify problem (initial values, parameter(s),...) Integrate further than 0.43368614

311

Restart Save results from 0.0000 to 0.43368614 Toggle out of quiet mode e x p e r t menu Quit [Our user selects ] B Enter up to I0 expressions to be plotted ..., (separated by commas and without embedded blanks) or <:RETURN> for your dependent variables against t, or ? for help, or * for preceding menu. x y [The user can plot any combinations of the variables and parameters. For this problem the phase plane (z,y) is needed, and the user responds to three further questions (about scales) with carriage returns, (omitted) to get the plot in Figure 3 below. On that we have also drawn in the tail section of the aircraft. ] [We mention that the user can plot any expressions using the variables or parameters against each other (up to nine simultaneously). For the case of our sample problem, we asked to plot z versus y, but could also plot any of the variables against t, perhaps also including rho*Cd*S*v*v/2, which is the air drag. Plots can be linear-linear, semi-log or log-log. In this as in almost all other cases, responding to a query vdth a carriage return will provide reasonable defaults ] [After the user has examined the plot on the screen he can obtain a copy on an attached printer by typing 'p'. PLOD also displays a cursor which can be moved about by using the keyboard's arrow keys and simultaneously reads out its position at the bottom of ~.h¢ ~crc¢n. This is useful for determining numerical values at points on the plotted curves. The user can also type 'z' which will draw a small box around the cursor and subsequently zoom into the box. PLOD's integration output points are used to drive its graphics. This differs from most other simulation packages in which the user must specify the plotting interval, or else the package chooses i~,s own equally spaced interval. Because the output points are unequal~,y spaced and chosen adaptively~ this avoids an aliasing problem that is common in generating plots of dynamically changing variables. Plots are guaranteed to be smooth and cannot miss structural details. We illustrate this in Figures 1, and 4, which are from the chemical oscillator problem mentioned earlier. To prevent too much data from being plotted, PLOD uses a "sieving" algorithm to only plot those points affecting the appearance of the graphs. Thus PLOD will use more of the original data points to generate Figure 4 than it uses for the equivalent region in Figure I. ] [After the plot our user is shown the Main menu again. This is the time to select M, for Modify problem, alter the value of parameter V¢, the initial conditions for v and O and finally integrate again. The data can also be saved on a file for other analysis.1 V L G r a p h i c s Capabilities We examined several Fortran callable libraries which were available for the IBM Personal Computer family, including IBM Graphical Kernal System [15], HALO [16], and GRAFMATIC [17]. Each had the capability and flexibility to satisfy our req'~rements. But the PLOD user needed to compile and link programs, hence the graphics library must be available to the linker. Distributing such a proprietary library would have violated each of the vendors' purchase agreements. Furthermore, we were committed to

312

development of public domain software. In this context graphics, or at least plotting, was seen as an important need. Consequently, we wrote a small set of screen plotting subroutines which were designed to be used with PLOD, and also could be called (from Fortran) by other scientists who were in a hurry to get their data displayed. These are also available to users who request them separately. VII. Conclusion One of the important results of the rapid spread of micro-computers is that people with less computer-oriented training are using computer.~ to help them in their jobs. This has pressured software developers into finding clever and interesting methods of generating problem solving tools which do not require programming. PLOD is one step in this direction. As our ingenuity and experience improve, we expect to see even better products with this same "non-programming" flavor. References 1. H. Fowler, D. Kahaner, J. Knapp-Cordes. F. Sullivan, "Wave Form Simulations for Josephson Juuction Circuits Used for Noise Thermometry" National Bureau Standarrl8 Teeh Journal 1983. 2. H. Kocak, Differential and Difference Equations t h r o u g h C o m p u t e r Experiments, Springer-Verlag, 1985. 3. Catalog of Simulation Software, in Simulation, Vol. 45, No. 4, pp196-209, i985. 4. G.R. Cook and E. Simiu, "Chaotic Motions of Forced and Coupled Galloping Oscillations," Sixth U.S. National Couference on Wind Engineering, University of Houstot, Houston, Texas, March 8 i0, 1989 5. D. Ka~aner. C. Moler, S. Nash, N u m e r i e ~ M e t h o d s and Software, Prentice Hall, 1989. 6. C. W. Gear, Numerical Initial Value P r o b l e m s in O r d i n a r y Differential Equations, Prentice Haft, 1971. 7~ ISSCO Graphics, 4186 Sorrento Valley Blvd, San Diego, CA 92121. 8. Hercules Computer Technology, Inc., 2550 Ninth Street, Berkeley, CA 94710. 9. Lahey Computer Systems, Inc., 31244 Palos Verdes Drive West, Suite 243, Rancho Palos Verdes, CA 90274. 10. Ryan-McFarland Corp., 609 Deep Valley Drive, RolLing Hills Estates, CA 90274. 11. Phoenix Computer Products, 320 Norwood Park South, Norwood, MA 02062. 12. R. C. Aiken, ed Stiff Computation, Oxford 1985. 13. Advanced Continuous Simulation Language (ACSL) - Reference Manual, Mitchell and Ganthier, Assoc., 290 Baker Av., Concord Mass. 01742, p. A-19, 1981. 14. D. Kahaner, W. Lawkins, S. Thompson, "On the Use of Rootfinding ODE Software for the Solution of a Common Problem in Nonlinear Dynamical Systems," Proc. Int. Congress on Comp. and Appl. Math., U. Leuven (Belgium), July 1988. 15. IBM Personal Computer Graphical Kerna] System, IBM Personal Computer Software, Number 6138807, Nov 1984. 16. HALO, Media Cybernetics Inc. 7042 C~rr,~l] Ave~ T~kcm~t Park MD 20912. 17. GRAFMATIC, Microcompatibles 11443 Oak Leaf Drive, Silver Spring MD 20901.

313

0 UO

J f

0 0

- -

1 I

Lf~

I

k. I

0

I I t I

u

0 0 04

I I

m

I

e m

I

u m o

I

o#_

I

Q

I

k_

I I

C7)

/

u

I

e m

I

u

I

I,

0 0

I

e B

I !

E ID

I

/

I

II

C)

/

I I

/

I I.t /

0 LO

/

I /

...

,

(2)

0 iJll

J l

liD

o .

,lIIIII

I

~

.

o

I

I,11111 i

~

.

o

.

,

I,lllli

l

I

UIIIII

I

i

i,,,,ll,

PO

04

~

o .

o

o

.

.

I

,mill

o

o

uo!loJlUeOUOO

I

i

illll,,

~

0

i

I

,,,,H,¶

o~

, 0

i

L~

~

, 0

[

314

0

i

0 m--

l

i

0 f~ r

0

A~

Om

0 ~a

0

o

°~

L_

0

o.--~

0

A.

X

W

o-

0 0

0

h-

m ou

0

0

i o ~0

(~

~-

0

O)

CO

~

~

~)

~

~0

(~

,-

0

315

Illlll

I

I

~llllll

I--

d,'rlll ,, I

Ii i

ulllll

i

i

1111111

g

i

ullll

i i

i

I I

C)

I I I I

Q. :3

I

cT~

00

0 m

I I

U

r'-.-

co

0

L

I f I

0 m

m

e m

u'3

I-\

m

\

0

c~

\ e ~

E

4) c-O

iJllll

O

i

i

llll¶lll

i

iJlllll

i

i

I, iiill

i

i

ll,llll

i

i

illlllll

i

o

U3

0

O

C~

~3 C7~ O m

I

r.O

E ew

0~

C)

Q)

0

T'-"

UO!IDJ~Ue~UO0

I