LIBRARY OF THE UNIVERSITY OF ILLINOIS AT URBANA-CHAMPAIGN srrv.o.S>S3-SS?> OOf> . Z CENTRAL CIRCULATION AND BOOKSTACKS The person borrowing this material is re- sponsible for its renewal or return before the Latest Date stamped below. You may be charged a minimum fee of $75.00 for each non-returned or lost item. Theft, mutilation, or defacement of library materials can be causes for student disciplinary action. All materials owned by the University of Illinois Library are the property of the State of Illinois and are protected by Article 16B of Illinois Criminal law and Procedure. TO RENEW, CALL (217) 333-8400. University of Illinois Library at Urbana-Champaign JO L ■*■ -i- LUU 1 JUL i I 2DD1 When renewing by phone, write new clue date below previous due date. L162 Digitized by the Internet Archive in 2013 http://archive.org/details/numericalsystems555brow YruXl ^ UIUCDCS- -R-T3- -555 Lz NUMERICAL ■ ROY SYSTEMS ON A MINICOMPUTER by LEONARD BROWN, JR. COO-1U69-0215 February 1973 UIUCDCS-R-T3-555 NUMERICAL SYSTEMS ON A MINICOMPUTER* by ROY LEONARD BROWN, JR. February 1973 DEPARTMENT OF COMPUTER SCIENCE UNIVERSITY OF ILLINOIS AT URBANA-CHAMPAIGN URBANA, ILLINOIS 6l801 Supported in part by the Atomic Energy Commission under contract US AEC AT(ll-l)lU69 and submitted in partial fulfillment of the requirements of the Graduate College for the degree of Master of Science in Computer Science. Ill ACKNOWLEDGMENT The assistance of the following people is greatly appreciated. My advisor, Professor C. W. Gear, provided both direction and cogent suggestions. Professors M. H. Pleck and R. L. Ruhl, and the Department of General Engineering provided the Illinois Graphics Computer System for my use in developing both the theory and the program presented here ; and my colleagues at IGCS, W. F. W. Tarn, D. Mueller, and T. Runge co-developed the software we all share . Mrs . Barbara Armstrong typed the manuscript. Finally, the Department of Computer Science and the U. S. Atomic Energy Commission supported the research. PREFACE This thesis defines the concept of a numerical system for a minicomputer and provides a description of the software and computer system configuration necessary to implement such a system. A procedure for creating a numerical system from a FORTRAN program is developed and an example is presented. The reader should have some knowledge of FORTRAN and minicomputer operating systems . A familiarity with PAL assembly language for the PDP-11 is necessary to make full use of the examples . V TABLE OF CONTENTS Page 1. THE PROBLEM 1 1.1 Definitions 1 1.2 Minicomputer Use 1 1.3 Outline 3 2. THE COMPUTER SYSTEM k 2.1 IGCS Hardware and Software 1+ 2 . 2 Alternate Systems 6 3. THE OVERLAY CONCEPT 8 3.1 Definitions 8 3.2 User-written Overlay Facilities 9 3 . 3 Manufacturer-supplied Overlay Facility 11 k. PARTITIONING THE FORTRAN PROGRAM lk k.l Criteria for Partitioning lk k.2 COMMON Variables lk k.3 Program Logic 15 k.k Control Transfer 17 5. USER SUPPLIED OVERLAYS 20 5.1 User Input 20 5.2 User Subroutines 20 5-3 Variable Size Arrays 22 6. CONCLUDING OBSERVATIONS 2k LIST OF REFERENCES 26 APPENDICES A. User -written Overlay Facility 27 B. DIF11 Program Code 30 1. THE PROBLEM 1.1 Definitions A numerical system as used here is any long program involving large amounts of arithmetic computation on a set of numerical data that is small enough that manipulation of the data is not a significant part of the problem. Manipulation means movement of the data from one storage device to another or change from one format to another, and similar programming problems not ordinarily associated with numerical analysis. Some numerical systems would be, for example, programs designed to solve difference equations, find roots of nonlinear equations, perform numerical quadrature of an input function, minimize a function of several variables, or, in the example dealt with here, to solve a system of ordinary differential equations. All of these programs have in common the require- ment of a great deal of arithmetic computation and simplicity of input: a number of polynomial coefficients, a difference equation of specified order and coefficients, a subroutine to evaluate a function with some initial conditions. 1.2 Minicomputer Use Such a numerical system is often difficult to write for a minicomputer for several reasons. Most minicomputers have very little computing power for the vast amount of arithmetic operations needed to run a program of this type; some do not even have a hardware multiply/ divide device, making complicated arithmetic operations time consuming. Small word sizes of 12 or 16 bits require that any floating point calculations use multi-word manipulations that take much time and storage space. The small size of main memory means that the program will probably not all fit into it at one time. This also limits the amount of numerical data since storing it anywhere hut main memory would be impossibly time consuming, given the block structure and slow access time of most peripheral storage devices. In spite of these drawbacks, there are still several reasons for creating a numerical system for a minicomputer. Although running a numerical system is time consuming, it is not as expensive as the same amount of time on a System/360 computer. Some minicomputers are idle part of the day and a set of numerical systems, possibly in a program library for the installation, would fill in the gaps by being run in free time or at scheduled times daily. Making such programs available could reduce an organization's expenses if these programs are normally taken to a large service computer. The complexity of multi-word manipulation, transcendental function evaluation, and programming of complicated formulae is reduced greatly because many l6-bit word size minicomputer software packages include both a library of subroutines that perform multi-word manipulations, and a FORTRAN IV compiler that allows execution of FORTRAN programs with approximately the same accuracy as a System/360 computer using single or double precision arithmetic. Finally, the problem of program size can be solved by proper use of peripheral storage devices such as a fixed head disc, drum, or other "random" access device. By storing large program segments that are always executed together on such a device (hereafter called disc although any "random" access device 'may be substituted), the program can be executed by bringing these segments into a core buffer, overlaying the program segment previously there. This core buffer will be called the overlay buffer. 1.3 Outline This paper will deal primarily "with a procedure for transforming a large FORTRAN numerical program known to work properly in a large computer environment into a numerical system that runs efficiently in a minicomputer environment. As a helpful example, a problem which has been solved using this procedure will be discussed wherever applicable. The coding of the example numerical system is found in Appendix B while the original problem is due to Gear [3] and involves integrating a first order system of ordinary differential equations through some interval of the independent variable, given a set of initial conditions and a subroutine to accomplish the differentiation. The original program, DIFSUB, is designed to handle stiff systems of equations as well as better conditioned ones. The topics to be presented are listed below. First, a computer system which is well adapted to the problem is described; criteria for choosing similar systems are listed as well as references to several such systems. Next, the concept of a program overlay is developed. The application of this concept in partitioning a FORTRAN program will be presented, and then several problems resulting from such a partition will be discussed. These last will fully describe the procedure and some additional comments about it will conclude the paper. 2. THE COMPUTER SYSTEM ; 2.1 IGCS Hardware and Software The Illinois Graphics Computing System used in the development of DIF11, the numerical system described in Section 1.3, belongs to the Department of General Engineering at the University of Illinois at Urbana-Champaign . As described in [5] and [6], it is based on a PDP-11/20 minicomputer from Digital Equipment Corporation (DEC). It has a byte addressable core memory of l6K l6-bit words with a ROM bootstrap loader. All communication between the CPU and peripheral devices takes place over a single bus which DEC calls UND3US. Access time for main memory is 490 ns. out and 8 ns . in. Instructions are basically of the two-address type, but 12 different addressing modes, eight general -purpose registers, and a software stack simulator allow a wide variety of instructions. Instruction formats require one to three words, depending on the addressing mode of operands. A KE11-A extended arithmetic unit provides single word hardware multiply and divide operations. A 256K word fixed head disc, accessible through the available software in 64 word blocks, has an average access time of 17 ms . and a transfer rate of 16 ys/word. A card reader, teletype, and Gould 4800 electrostatic printer /plotter complete the basic system. The Gould 4800 is used by DIF11 as a printer through a software program, LP, since the Gould will eventually be the output device of a system containing both a CalComp-type graphical language and a simulation system of which DIF11 will be part [5 ] • The software package available with a minicomputer is very important since it can greatly reduce the work needed in programming a large system. IGCS currently operates with a Disc Operating System, DOS version UA, which includes utilities for transferring data "between memory and peripheral devices, an assembler, a FORTRAN IV compiler, and a relocating Linkage Editor which presently includes an overlay facility. See [2] for details. EAELIB is a library of FORTRAN functions (e.g. SIN, COS) and utility routines (e.g. floating point addition) designed to operate with the KE11-A multiply /divide unit. These routines are linked with a FORTRAN program by the Linkage Editor and are used by the program at execution time. This system is suitable for a numerical system since l6K of core is adequate for keeping the operating system core (DOS is not entirely resident, but overlays part of itself depending on its current operation), a set of routines from EAELIB for a large program, a main program, the overlay buffer(s), and storage for data whose size is variable and dependent on an input parameter. The FORTRAN compiler is very efficient in compiling code to do the actual computation involved, although control transfer, allocation of memory space, etc., are best performed using assembly language interfaces. FORTRAN programs interpret a polish string of addresses of routines placed in core by the Linkage Editor from EAELIB; these strings are the result of compilation of the program. Storage is best handled as FORTRAN arrays whose elements fill consecutive four byte segments of core and are easily processed by either FORTRAN or the assembly language interfaces. The original working version of DIF11 was implemented in 12K core with user-written overlay software so the present system is more than adequate. The speed of access, relatively small blocksize, and adequate device driver software make the fixed head disc easily accessible for transferring unchanging program segments to core, 2.2 Alternate Systems - Inspection of the above system leads to a set of criteria which a minicomputer system should satisfy to create a well suited environment for a numerical system. These criteria follow. The main memory should he an adequate size to contain the system's resident monitor, a large section of code from the numerical system, and all the variable storage for the system. An adequate word size to minimize multi- word floating point manipulation is important; at least 16 bits is recommended. The execution speed should be fast enough to allow an adquate turnaround time for the size problem being run; a good criterion would be that the minicomputer turnaround time including loading the numerical system and execution should be no more than 10 times the turnaround for a full size service computer from input to retrieval of output. The peripheral device that acts as a bulk memory, usually a fixed head disc or drum, should have fast transfer into core. Any sequential device would be unsuitable for this . An efficient FORTRAN compiler (i.e., programs compiled use close to the absolute minimum of core) is necessary. Since many numerical systems may need to use double precision arithmetic, this feature should compile and execute efficiently. Finally, a good software monitor capable of fast and effective control of device transfers while occupying a minimum amount of main memory is important, both at execution time and during creation of the numerical system. Several computer systems, while dissimilar to the PDP-11/20 in many ways, appear to meet the above criteria. Some may be better suited to the purposes of the purchaser than others and still provide an acceptable environment for numerical systems. They are: Hewlett-Packard 2100 A [k] Interdata Model ?0 [7] Systems 72 [8] Varian 620/f [9] Westinghouse 2500 [10] 3. THE OVERLAY CONCEPT 3.1 Definitions A program overlay is a segment of executable code intended to be used with a permanently resident program segment. However, this code is stored permanently only on a peripheral storage device (disc) and is brought into a core buffer only when needed. When another overlay is needed instead, that overlay is brought from disc and placed into the same buffer, thus "overlaying" the previous segment. Clearly, N overlays of approximately equal size effectively expand the overlay buffer to simulate a memory area N times as large. However, one pays for this in terms of the cost of disc to core transfer time. This is more advantageous than having all executable code in core (if possible) and using the disc to store variables because these must be transferred not only into core, but also back from core to disc after being processed, whereas overlays are constantly "refreshed" from an unchanging disc image. Since overlays for numerical systems are usually parts of a FORTRAN program originally intended to be executed in the same main memory, the overlay should have access to all variables stored in resident core, and to all resident subroutines, especially such utility routines as are normally used by minicomputer FORTRAN compilers for complicated arithmetic operations such as floating point addition (e.g. EAELIB in Section 2.1). Because of these features , certain control paths between the resident code and each overlay must be provided. These paths are: a. Resident to overlay b . Overlay to overlay c. Overlay to resident; return of control d. Overlay using resident subroutines A core map of an idealized numerical system is in Figure 3.1.1. MAIN MEMORY DISC oiuii FLOATING POINT UTILITY ROUTINES MAIN RESIDENT FORTRAN ROUTINES OVERLAY BUFFER VARIABLY DIMENSIONED ARRAY STORAGE (STACK) ~f~ 000000 DOS MONITOR FLOATING POINT UTILITY ROUTINES OVERLAY 1 OVERLAY 2 Figure 3.1.1. Sample Program Map with Control Paths Indicated 3.2 User-written Overlay Facilities Some minicomputer software systems do not have overlay facilities; in these cases the user must code his own. Such a system was written for use with IGCS before versions of FORTRAN and the Linkage Editor were provided with one. The code for this system in PAL assembly language is in Appendix A, and a core image for creation of an overlay file with it is in Figure 3.2.1. 10 r FLOATING POINT UTILITY ROUTINES FINAL J LINK \ MAIN-FORTRAN ROUTINES CONCTM OVERLAY MAKOVL (STACK) "J"" 000000 DOS MONITOR } OVERLAY BUFFER DISC Figure 3.2.1. Core Map for Creating an Overlay Using MAKOVL The system consists of two programs called CONCTM and MAKOVL. CONCTM is a core resident subroutine which contains tables of the file names of all overlays, their length in words, and their starting addresses in core. On being called, it looks up the disc start block of the overlay using the FILE BLOCK and LINK BLOCK [2] and places these in the TRAN BLOCK. The overlays are contiguous on disc, i.e., if they occupy more than one disc block, they are placed in consecutive disc blocks. When an overlay is called, the call is to CONCTM with the first argument the integer designating the overlay wanted and the following arguments being the proper arguments for the subroutine from which the overlay was made. CONCTM transfers the overlay into core at the indicated starting address, using a TRAN request and the TRAN BLOCK mentioned above. Since the addresses of the arguments of the overlay subroutine follow the statement JSR CONCTM, R5 in the main program, and since subroutine arguments are referenced relative to the return address stored in register R5, a simple transfer to the start 11 address of the overlay starts its execution. A FORTRAN RETURN statement shifts control back to the resident program segment. The overlays are created "by linking together all resident parts of the program, the overlay, and MAKOVL below these (the Linkage Editor links in the order of its input string, filling the top of core). MAKOVL is executed from the keyboard and computes the length in words and disc blocks of the overlay and transfers the overlay to disc under the name placed in its FILE BLOCK (this can be set from the keyboard) . The user saves the length and starting address and inserts them into the final version of CONCTM when the routines indicated as the final link in Figure 3.2.1 are linked together as the initial load module. Here it is difficult to separate the floating point package routines between those needed by the resident segment and those needed only by the overlays . One solution is to determine all utility routines needed by any program and force them to be permanently resident. On the PDP -11/20, this can be done by linking each program and inspecting its load map which lists all utility routines called by it. The union of all utility routines listed on any load map is then placed in a GLOBL state- ment. A version of DIF11 using the above system was programmed and tested. 3.3 Manufacturer-supplied Overlay Facility The overlay facility described here is found in more detail in [l], It is typical of those available from manufacturers' software groups. The system uses a single multiply entrant subroutine LINK with entry points LINK and RETURN. These provide all the control paths described in Section 3.1. CALL LINK('FILE') in resident FORTRAN code or JSR LINK,R5 with 12 the proper argument in assembly language l) initializes FORTRAN traceback routines, 2) saves register 0-5, 3) causes the named file to be found, transferred into core, and executed as if it were a main program. The statement CALL RETURN restores register 0-5 and returns control to the place in resident code where LINK was last called. Since it is easiest to write and test the FORTRAN part of the numerical system with the overlays as subroutines, a problem arises in that the overlays must be set up as main program segments . In some cases the FORTRAN SUBROUTINE statement can be removed, but if any parameters must be passed as if the overlay were a subroutine (see Section 5«2), then an assembly language main program can be written with the appropriate arguments. Such an example is the routine INTRO in Appendix B. 011111 INITIAL LINK 000000 r RESIDENT FLOATING POINT UTILITY ROUTINES V MAIN -FORTRAN ROUTINES INTRO OVERLAY FLOATING POINT UTILITY ROUTINES OVERLAY SUBROUTINE (STACK) T DOS MONITOR "S OVERLAY BUFFER DISC Figure 3.3.1. Core Map for Creating an Overlay Using Manufacturer-supplied Software The overlays are created using the relocating Linkage Editor, LINK. The resident program segments are first linked together, bringing 13 all floating point utilities used by the resident code into core. Two of the outputs of LINK are saved: the object module, and the Symbol Table -which lists all subroutine start addresses including those for the utility routines. Then each overlay (with its main program INTRO if needed) is created; the first input to LINK is the Symbol Table, and the top address of the output load module is set to two less than the bottom of the core resident module just created. The Symbol Table input allows the overlay to use utility routines already core resident and only utility routines which it alone needs are added to the overlay. Figure 3.3.1 shows the configurations of these different uses of LINK to create the resident and overlay modules . Ill h. PARTITIONING THE FORTRAN PROGRAM k.l Criteria for Partitioning After providing for an overlay facility, one next considers the overlays themselves. It is best to begin -with a FORTRAN program which is known to run in a large computer environment and insure that the overlay interconnections will work by breaking the program down into overlay subroutines and running it again on the same computer. There are really no hard and fast rules that can be applied, but since we are dealing with relatively slow minicomputers and one of the slower operations is disc-to-core transfers, one goal of the partitioning process is to minimize the number of such transfers . Since core space is at a premium, most of the original program should be put into overlays, leaving minor bookkeeping arithmetic and major control transfer decisions in the resident program segment; by making the overlays of nearly uniform length, the overlay buffer (which is as long as the longest overlay) is also made as small as possible. The implementations of these two goals are clearly at odds since the real minimum number of overlay transfers from disk would maximize the overlay buffer, and the minimum overlay size would be only two or three FORTRAN statements, requiring innumerable disc transfers, usually of length less than one disc block. Proper implementation calls for a balance between these two extremes. U.2 COMMON Variables Because the contents of the overlay buffer at the time a new overlay is transferred in will be overwritten, any FORTRAN variables 15 (other than special constants) must be core resident. This is most easily accomplished with COMMON statements which, when linked into load modules by most Linkage Editor facilities, leave all of the variables in a resident core block. When the overlay modules are linked, the addresses of elements in the COMMON block are available from the Symbol Table and are linked with the overlay. To form the COMMON block, each overlay and the resident program is inspected and all of its fixed length variables are listed. A variable array is of fixed length if its dimension does not depend on an input to the program. Next, any variable appearing on any two lists must obviously be shared by two different routines and is placed in the COMMON block. Any data that remains the same in an overlay and is not shared by any other overlay can be defined by a DATA statement and need not be in the COMMON block. It may be convenient to have several COMMON blocks since some variables may be used exclusively by certain overlays and no others . Allocation of space for arrays with variable dimensions will be discussed in Section 5.3, but note that the initial address of each such variable must be passed to the overlay as a subroutine parameter (it cannot be in a COMMON block since its initial address is not available at Linkage time), and the variable dimension must be passed to the overlay as a subroutine parameter as required by most FORTRAN compilers . h.3 Program Logic Some few FORTRAN programs are executed sequentially, i.e., starting at the first statement and "dropping through" to the end with a few tight DO-loops . Partitioning these is trivial. The program is divided into equal 16 length overlay segments; either each overlay calls the next, or, if this is not possible, a short resident program calls each overlay in sequence. However, most programs are not this easy. All or part of the more usual program is executed iteratively either a fixed number of times or until certain conditions are met. DIF11, for example, predicts next values of the differential variables, corrects these values, and, under certain conditions, changes other control variables (order, step size, etc.). This whole cycle — predict, correct, change — is repeated until a final value is reached for the independent variable. The best approach to partitioning such a program is to follow a principle of least use; if something is seldom used, it is an ideal candidate for an overlay. Many programs have a long initialization section which is only executed the first time the program is used. Clearly, this should be an overlay. The original version of DIF11 had such an overlay but the present version does not need one since the initialization was shortened. A segment that is executed at most once per iteration is the logical next candidate. In DIF11, SUB3 and SUB5 are executed only once per iteration and SUB1, SUB2, and SUB6 are executed only when a change of control variables is being considered or carried out. (See Appendix B). DO-loops are best placed completely within an overlay to avoid each pass through the loop causing one or more disc transfers to the Overlay buffer. Occasionally, however, a DO-loop will contain enough material for several overlays (remembering that uniform overlay size is a goal of the partition), and the DO-loop then contains several overlays. IT DIFFUN, SUBT, and SUBU are all contained in one loop in DIF11. The same loop contains code for the numerical evaluation of the Jacobian matrix which could have been an overlay. It was not because this would have required one overlay to call another and to return to the original overlay. This complicated process is best avoided but could be handled by using two overlay buffers. See Section 5.2. Occasionally, one especially long overlay may be needed but only part of it need be in core the whole time it is being executed. In such a case only enough of the overlay should be brought into core at one time to fill the regular sized overlay buffer, and when some part of it that will not be executed again is completed, the rest is overlayed there. This will be possible if the overlay facility allows overlays to call overlays as specified in Section 3.1. Any segment that is called from a number of different places in the resident code is probably best left in core, but if it is long or called very seldom, it can be made an overlay as in SUB8 in DIF11. k.h Control Transfer It must be noted that each overlay may have several GO TO statements which reference a statement in the resident module, or worse, in another overlay. The partition should minimize this, but if it happens, than a control transfer vector can be used to allow this. Either a single variable, JLINK, or an array, JLNK(N), can be used; the various values should have some significance to the program. The single variable is probably better since it is less complicated. Set JLNK = 1 on entering 18 an overlay, then if a statement is, say GO TO 750 and 750 is in the resident module, change the statement to JLNK = 2 RETURN If a statement is IF( LOGICAL EXPRESSION) GO TO 715 change this to IF (LOGICAL EXPRESSION )JLNK= 3 IF ( JLNK. EQ. 3) RETURN This requires that any arithmetic IF statements he made logical IF statements Immediately after the overlay calling statement in the resident code, a computed GO TO statement can he used to effect the transfer: CALL OVERLAY (SUB 5) GO TO(5 1 40,750,715,.. .),JLNK 5^0 ... remembering the possibility of executing the next sequential statement if JLNK remains 1. Transfers into the middle of another overlay can be handled by a computed GO TO at the beginning of the overlay so that this is the first statement other than subroutine calls after the computed GO TO in the resident code. JLNK must be set to 1 after executing the initial computed GO TO statement. For example, GO TO TOO 19 "becomes JLNK=2 RETURN * * * _ ********** 5U0 GO TO(550,600),JLNK 600 CALL OVERLAY (SUB2) *** ********** ( overlay 1 ) (MAIN) GO TO(650,TOO) ,JLNK TOO JLNK=1 (overlay 2) Examples of all these techniques are in SUB5 of DIF11 in Appendix B 20 5. USER SUPPLIED OVERLAYS 5.1 Us er Input Many numerical systems can operate with a minimum of input. For example, the only input needed to solve a difference equation is the order, the coefficients, and the initial values. These can all be input on cards or keyboard, and the storage for this data can be in the core resident main program by specifying a maximum size for each variable array. Other programs will require more from the user than numerical data; a numerical quadrature program will require that some way be provided to evaluate the function at any given point, usually a subroutine Some programs will have large data arrays which differ greatly in size depending on input parameters . Storing these within the main program segment may be too restrictive. An LU matrix decomposition would need such a variable sized space to store the matrix and outputs. And, some systems require both an evaluating subroutine and variable sized arrays, as is the case with DIF11. 5.2 User Subroutines The solutions to these problems for DIF11 are presented here; they seem applicable to similar cases on other systems and are at least a good starting point. In DIF11, the user is responsible for writing a FORTPAN subroutine called DIFFUN(T, Y,DY) , and by following carefully detailed instructions, creating an overlay from it. First, compile DIFFUW(T,Y,DY) where T is the independent variable, DY(l) is the first derivative of Y(l), and all derivatives of order m higher than one have been replaced by m first-order equations using standard techniques. 21 Y(l) and DY(l) should be dimensioned N, the number of equations. Next, call LINK, the Linkage Editor. Specify DIFFUN. LDA as output; specify as input ST (the Symbol Table of Section 3.3 which is kept available on disc or tape), the object module of INTRO (also available), the object module resulting from compiling DIFFUN, and the EAELIB. Also, specify that the top core address available to this module is two less than the bottom of the resident code (this address is posted in the machine room). After running LINK, run the program DIF11 with two data cards: one specifies N, T , H (the length of the integration), EPS (the error criterion), and ISTORE (the octal top address to be filled by dynamically allocated variable arrays); the other lists the N initial values of Y(l). This overlay uses the same overlay buffer as all other overlays . This is possible because none of the other overlays call it. Any other overlay could call DIFFUN and then return to the main program; however, a great deal more work would be required to let an overlay call DIFFUN, have DIFFUN reinstate that overlay, and then transfer to the statement after the call to DIFFUN. This would be equivalent to a regular subroutine call to DIFFUN. Another possibility is to provide a second buffer just for the user supplied overlay, allowing all overlays to call it. However, it must be remembered that this takes valuable core that might otherwise be used for variable storage. Since there is no way to predict the length of a user supplied module, this allows the first overlay buffer to be of fixed length while the second varies. This may be an important consideration for some users. 22 5.3 Variable Size Arrays The problem of variables whose dimensions depend on an input parameter is best solved by allooating storage space for the FORTRAN arrays in an assembly language MAIN program. In DIF11, various arrays are dimensioned (N,T), (N,2), or (N,N) . The MAIN program (see Appendix B) computes U*N (each floating point -word uses four bytes on the PDP-11/20), 28*N, 8*N, and U*N**2. Using RO (register 0) initially set to ISTORE as a pointer, MAIN subtracts an appropriate number of bytes from RO and assigns RO to an address vector ADDR . Then it places the addresses of the appropriate variables in the address list of the calling sequences of INPUT, OUTPUT, and DIFSUB— the FORTRAN subroutines called from MAIN — and executes. It should be noted that most FORTRAN compilers require a variably dimensioned array, e.g. Y(N,7), to be in a subroutine with the variable dimension N as an argument, so DIFSUB, INPUT, and the overlays must be subroutines. Whenever DIFSUB calls one of its overlays , the overlay starts with a main program INTRO as required by the overlay facility, and INTRO places the initial addresses assigned to the variably dimensioned arrays and kept in an assembly language named .CSECT block (same as FORTRAN COMMON) into the FORTRAN calling sequence. Then INTRO executes the subroutine call, returning to the overlay main program before calling RETURN (see Section 3.3) in order to restore all registers. In this way, eight variably dimensioned arrays are allocated dynamically in a manner similar to ALGOL by using assembly language interfaces between FORTRAN program segments . It should be noted that the pseudo-stack mentioned in Section 2.1 is kept immediately below the overlay buffer by the DEC provided software, 23 and extends into lower core addresses . This requires that the dynamically- allocated variables he placed well below the stack and is the reason ISTORE is an input parameter. If a user written overlay facility were available, the stack could be moved at will below the data. 2k 6. CONCLUDING OBSERVATIONS The following observations seem appropriate to the problem discussed and are presented in an attempt to provide some perspective. Firstly, the method of reducing a large FORTRAN program to a mini- computer overlay system is best considered as a procedure. A list of step by step instructions can be found in the preceding chapters , but no guarantee is provided that the procedure will ever stop. Thus, this is not an algorithm with input: one large FORTRAN program; output: equivalent minicomputer numerical system. A counter example is a program with larger storage requirements than could ever fit into core or be profitably allocated to disc. This procedure can be helpful but no promises are made . Secondly, in most cases the more general a program is, the larger it becomes. A simple Runge-Kutta (R-K) program with step doubling will take much longer to integrate a complicated system on a 360/75 than DIFSUB will; but if a problem is fairly simple, a R-K program requiring no overlays would be much faster than DIF11 on the PDP -11/20, Thus, one can save time by having available a simple method for simple problems as well as a complicated method to solve complicated problems. Finally, one should be wary of manufacturer supplied overlay facilities, especially if overlays are called more than a few times. The overlay processor supplied by DEC looks up the current disc start block of each overlay every time it is called, even though some overlays are called hundreds of times' from the same disc location. For this reason, the final simulation language package, ILLISM, which DIF11 will be part of, will use a locally written overlay facility that looks up disc addresses once, stores them in a directory, and thus executes twice as fast. 25 With these comments in mind, one can conclude that many large numerical analysis programs can he profitably used on a minicomputer system using the techniques described here. 26 LIST OF REFERENCES [l] Digital Equipment Corporation, "Getting FORTRAN on the Air," (DEC-11-SFDC-D) , Maynard, Massachusetts, 1972. [2] Digital Equipment Corporation, "DOS Monitor Programs Handbook," (DEC-11-SERA-D) , Maynard, Massachusetts, 1971. [3] Gear, C. ¥. NUMERICAL INITIAL VALUE PROBLEMS IN ORDINARY DIFFERENTIAL EQUATIONS, Prentice-Hall, Inc., Englewood Cliffs, New Jersey, 1971. [k] Hewlett-Packard Company, "2100A Computer Reference Manual," (102100-90001), Cupertino, California, 1970. [5] Pleck, M. and Ruhl, R., "Illinois Graphics Computing System," Progress Report #1, Department of General Engineering, University of Illinois at Urb ana-Champaign , 1971. [6] , Progress Report #2. [7] Interdata, Inc., "Model 70 User's Manual," (29-26l), Oceanport, New Jersey, 1971. [8] Multidata, Inc., "Computer Reference Manual — Systems 72," (31101), Westminster, California, 1970. [9] Varian Data Machines, "Varian 620/f Computer Handbook," (98-A 9908 001), Irvine, California, 1970. [10] Westinghouse Computer and Instrument Division, "Computer Reference Manual," (25 KEF -001) , Orlando, Florida, 1971. 27 APPENDIX A User-written Overlay Facility 28 QNO NUMBER OF WORDS. PSR NUMSEQ fiSR NUMSEQ ODD #1 , NUMSEQ ALLOC LNKBLK, FILBLK.NUMSEG MOV NUMSEQ, - MOV #FILBLK,- MOV #LNKBLK, ~ EMT 1 5 INC + BN E ERR .OPENC LNKBLK, F1L9I.K JNE ERR .OPENC LNKBLK, Fll MOV #13,FILBLK-2 MOV #FILBLK, - MOV #LNKBLK, - EMT 16 MOV #BUF,- MOV #LNKBLK,- EMT 2 WRIT LNKBLK MOV #LNKBLK, - E" M T 1 .CLOSE LNKBLK MOV *LNKBLK, -(8P> EMT 1 7 ERR: OVRLEN: .WOKD NUMSEQ: .WORD L N N K BLOCK •WORD .RLSE LNKBLK MOV #LNKBLK, - EMT 7 .EXIT EMT 60 . WORD .WORD 0,0,. • DF/ 1 LNKBLK: .RPD50 ' vr r ; FILE BLOCK . WORD 0, F1I.BLK: .RRD5B 'DIF' .RRD50 •FUN/' .RAD50 'LDA/ .BVTE 0,0,0,0 BUF: .WORD OVRLEN, 7, OVRLEN, I . END 'IFFUN 29 SUBROUTINE CONCTM •TITLE CONCTM SP « 16 CONCTM: POINTR: SP « *6 R0 " f e GLOBL CONCTM S2<*5> .POINTR POINTR POINTR R0, - MOV OSL ASL MOV MOV MOV MOV MOV MOV MOV MOV EMT MOV CLR MOV EMT MOV TST TST MOV MOV EMT MOV EMT MOV EMT JMP POINTR, R8 5 PUT 9EC0I POINTR , FILBLK+2 JFJLE Bl P0INTR-2, TILBLK ; PUT ADDR , TRNBLK*4 ADDR-2CR0) , TRNBLK+2 {PUT +,R0 JMOVE Sft jPUT 4*POINTR WORD OF ON INDICATOR TO POINTR THE STACK. MODULE 5045 IN FILE BLOCK TRAN BLOCK SECOND " BLOCK FIR9T WORD ( SP) ♦ , R0 iLNKBLK, - , TRNBLK STORT LOCATION IN •FILBLK - * 5 e . RAD50 . R ft D 5 O . R0D50 . RAD50 . RAD50 . ROD50 . R H D 5 * R ft "> NUMWRD:ADDR: . WORD . WORD LNKBLK: .WORD . RRD30 . WORD F ILBLK: . WORD . RAD50 . B V T E TRNBLK: .WORD . END #TRN8LK, - #LNKBLK, -<8P> 1 ttLNKBLK, - 1 *LNKBLK, -< SP) t STRNBLK+2 ;PO y D I F / /FUN/ /PAS/ /COL/ /QFI/ /ND/ /SCfi/ /LE/ /EVfi/ /L/ /CHB/ /NGE ' /INI/ /T/ . I N IT ; PUT .WRITE . WRIT . ALSE • LOOK START BLOCK IN TRAN BLOCK 5 .RLSE J 00 TO SUBROUTINE STORT INTER BND LIST OF FILE NAMES ADDRESS /COR/ /ECT/ WORD 0, 37260, 76,37122, 155,371 12, 161 ,36476, 36634,318,3 9064, 1 174,37072, 161 ,36414,420 ;link block 0,0, 1 /DF/ 0,4 ; FILE BLOCK 367 0,0. 1 /DF/ 0,4 0, /LDA/ 0,0,0,0 0,0,0,4 30 APPENDIX B DIF11 Program Code 31 THIS PROGRAM IS THE MRIN PROGRAM OF 1 F t 1 . IT INITIALIZES CONSTANTS, READS IN VARIABLES, AND ASSIGNS STORAGE LOCATIONS TO VARIABLY DIMENSIONED STORAGE FOR DIFSUB. .TITLE MAIN Re - se R1 R2 R3 R4 R5 SP PC f 1 *2 I 3 *4 *5 X6 .GLOBL . CSECT BEGIN: LOW), L0W2, OUT PUT, INPUT, LINK, RE TURN, N, DIFSUB .WORD a RET1 : ON: NORD LISTO: RET2: LOOPS: L00P4 ST 1 9: ST20: LIST: out: OUTLIS: RET3: N , INPUT RET1 , L.0W2, XLIM LOW1 , L0W2 ON L0W2, LOW I N, R1 R1 R1 LOH1 , R0 • 2, RO R I , R0 R8, ADV Rt , R0 RB.AERROR UB R 1 , R0 R0.AVMAX R1 , R0 RB, AI R1 , R3 R3 R3 R3 R1 , R3 R3, R0 R8, AV R3, R0 R0, ASQVE R1 R1 , R0 R0, ACSAVE N2, R0 R9, APM ♦ 4 PIP CLR SR R5 BR ,N,N2 CMP BMI MOV MOV ASL ASL MOV SUB SUB MOV SUB MOV S MOV SUB MOV MOV ASL ASL ASL SUB SUB MOV SUB MOV ASL SUB MOV SUB MOV BR MOV AV, LISTO JSR R3, INPUT ER RET2 . WORD 0, N, 0, 0, MOV AVMAX.R1 MOV N,R0 MOV #40200, CR1 )♦ CLR ♦ CLR ( R1 SUB # 1 , R0 CMP #0,R0 BNE LOOPS MOV AV, OUTLIS MOV #20, R3 MOV #ADDR,R0 MOV #LIST,R1 MOV (RB) + , (R1 >■» TCT CDOS- L U U K MOV AV.OUILI MOV #20, R3 10V #ADDR,R0 .10 V • L I S T " MOV * TST - C M P # , R 3 BNE L00P4 CLR JSTART BR ST20 MOV #1 , JSTART JSR R5, DIFSUB BR OUT .WORD N . WORD 0, 9, *, 8, 0, e, JSR R3, OUTPUT BR RET3 . WORD 0, N, XLIM CMP N t , KFLAG BEQ 3T1 9 r u t /.a ;CALL INPUT SUBROUTINE WITH N-0 TO READ IN JN, T.H.EPS, AND L0W2-BOTTOM ADDRESS OF {DIFFUN OVERLAY PROVIDED BV USER. [LOWEST OVERLAY START ADDRESS GOES TO L0W1 JTO BECOME TOP OF STORAGE AREA. {COMPUTE ADDRESS OF DYNAMICALLY ALLOCATED {STORAGE FOR DIMENSIONED VARIABLES. ;R1 CONTAINS 4*N, R0CONTAINS TOP OF STORAGE. JADV - • OF DY ;AERROR - # OF ERROR JAVMAX «■ • OF VMAX(N> j AIP • • OF IP ', R3 CONTAINS 7 + 4*N ; AV - # OF Y ;ASRVE • • OF SAVE(N,7) ;R1 CONTAINS 0*N {ACSAVE • # OF CSAVE*4 JAPW - « OF PW(N,N) ; READ V(*,I) FROM CARDS ; SET VMAX<*> - 1.0 JR0 18 THE LOOP COUNTER {PREPARE PARAMETER LISTS FOR CALLING {DIFSUB AND OUTPUT { JSTART { JSTART 6,0 32 N - .word e,0 N2: • WORD 8, 8 LOM1 : . WORD 03351 2,0 L0W2: .WORD 0,0 XLIM: .WORD 0,0 .CSECT COM JSTPRT - .+40 KFLOG ■ .+50 T = .+74 H = .+124 HMHX - .+140 HMIN » .+144 .CSECT flDDR HC'DR:flDV: .WORD PERROR: .WORD flVMflX: .WORD BIP: . WORD a v : . w R D Q pspve: .word pcspve: .word PPW: .WORD .END BEGIN Resident I/O Routines .T.TOLD, OWN, ENQ1 ■ FlflQ.NQ.NQOLD.NEWQ, EPS, T, TOLD, HNEW,HMQX,HMIN,E,EUP,EDWN,ENG>1, Resident FORTRAN Main Program 33 r SUBROUTINE DIFSUB COMMON /COM/- A<7>,JLNK, JSTART, K,KFLAG,NQ,NQOLD,NEWC!,EPS,T, TOLD, 1 IDOUB, IWEVAL, I RET , IRET1 ,H, HOLD, HNEW.HMAX, HMIN, E.EUP.EDWN, EN C'1 , 2ENQ2,ENQ3,BND,BR,DEL,DEL1 , D r + *»* + + + + + * + + * + + + + + *** + ** + + **** + * + * + + * + * + * + *4-*** + + + + ** + *±** + + + ** + + .t. + + 4-** c * * C* THIS SUBROUTINE INTEGRATES A SET OF N ORDINARY DIFFERENTIAL FIRST * C* ORDER EQUATIONS OVER ONE STEP OF LENGTH H AT EACH CALL. H MAV BE + C* INCREASED OR DECREASED WITHIN THE RANGE HMIN TO HMAX TO C* ACHIEVE AS LARGE A STEP AS POSSIBLE WHILE NOT COMMITTING A SINGLE + C + STEP ERROR WHICH IS LARGER THAN EPS IN THE L-2 NORM, WHERE EACH + C* COMPONENT OF THE ERROR IS DIVIDED BV THE COMPONENTS OF VMAX. + C + c* THE PROGRAM REQUIRES THE SUBROUTINES NAMED * C*DIFFUN * C* DECOMP(N,M,PW,IP) * C * SOLVE'N H PW CSAVEd I) IP) * C* THE FIRST, 'dIFFUN, EVALUATES THE DERIVATIVES OF THE DEPENDENT + C* VARIABLES STORED IN V < 1 , I > FOR I - 1 TO N, AND STORES THE * C+ DERIVATIVES IN THE ARRAV DV. DECOMP IS A * C* STANDARD LU DECOMPOSITION WITH PIVOTING THAT DECOMPOSES THE MATRIX * C* PW, LEAVING THE PIVOTS IN THE INTEGER ARRAV IP. M IS THE DECLARED * C* SIZE OF PW. IP IS SET TO IF PW IS SINGULAR. SOLVE PERFORMS * C * BACK SUBSTITUTION ON THE CONTENTS OF C S A V E < I , 1 > , LEAVING THE * C* RESULTS THERE. + C * * C* THE TEMPORARV STORAGE SPACE IS PROVIDED BV *MAIN* IN THE * C+ INTEGER ARRAV IP, THE SINGLE PRECISION ARRAVS PW, DV, + C* PRECISION ARRAVS SAVE AND C8AVE. THE ARRAV PW IS USED ONLV TO HOLD + C* THE MATRIX OF THE SAME NAME, AND SAVE IS USED TO SAVE THE VALUES * C* OF V IN CASE A STEP HAS TO BE REPEATED, BUT CSAVE IS USED TO HOLD * C* SEVERAL ARRAVS. * C + * C* THE PARAMETERS TO THE SUBROUTINE DIFSUB HAVE * O THE FOLLOWING MEANINGS.. + C+ * C* N THE NUMBER OF FIRST ORDER DIFFERENTIAL EQUATIONS. N + C* MAV BE DECREASED ON LATER CALLS IF THE NUMBER OF * C * ACTIVE EQUATIONS REDUCES, BUT IT MUST NOT BE * C* INCREASED WITHOUT CALLING WITH JSTART -0. * C* T THEINDEPENDENTVARIABLE. * C* V A 7 BV N ARRAV CONTAINING THE DEPENDENT VARIABLES AND * C * THEIR SCALED DERIVATIVES. Y < J ♦ 1 , I > CONTAINS * C* THE J-TH DERIVATIVE OF V SCALED BV * C * H + + j'FACTGRIAL<'J> WHERE HIS THE CURRENT * C * STEP SIZE. ONLV V<1,I> NEED BE PROVIDED BV * C + THE CALLING PROGRAM ON THE FIRST ENTRV. * C * H THE STEP SIZE TO BE ATTEMPTED ON THE NEXT STEP. + C* H MAV BE AD. JUSTED UP OR DOWN BV THE PRO GRAM * C* IN ORDER TO ACHEIVE AN ECONOMICAL INTEGRATION. * C* HOWEVER, IF THE H PROVIDED BV ♦MAIN* DOES * C* NOT CAUSE A LARGER ERROR THAN REQUESTED, IT * '"* WILL BE USED. TO SAVE COMPUTER TIME, * M A I N * * C* USESPFAIRLVSMALLSTEPFORTHEFIRST * C * CALL. IT WILL BE AUTOMATICALLY INCREASED LATER. * C* HMIN THE MINIMUM STEP SI2E THAT WILL BE USED. * C* UMAX THE MAXIMUM SIZE TO WHICH THE STEP WILL BE INCREASED * C* EPS THE ERROR TEST CONSTANT. SINGLE STEP ERROR ESTIMATES * C * D I V I D E D B V Y M A X < I > MUSTBELESSTHANTHIS * C * IN THE EUCLIDEAN NORM. THE STEP AND'OR ORDER IS * C* ADJUSTED TO ACHEIVE THIS. * C * VMAX AN ARRAV OF N LOCATIONS WHICH CONTAINS THE MAXIMUM * C* OF EACH V SEEN SO FAR. IT SHOULD NORMALLY BE SET TO * C * 1 IN EACH COMPONENT BEFORE THE FIRST ENTRV. 'SEETHE * C* DESCRIPTION Or EPb. ^ + C * ERROR AN ARRAV OF N ELEMENTS WHICH CONTAINS THE ESTIMATED * C* ONE STEP ERROR IN EACH COMPONENT. * C * KFLAG A COMPLETION CODE WITH THE FOLLOWING MEANINGS.. * C* +1 THE STEP WAS SUCCESFUL. * C* -1 UNRECOVERABLE ERROR + C* INTEGRATION COMPLETED * C* JSTART AN INPUT INDICATOR WITH THE FOLLOWING MEANINGS.. + C* PERFORM THE FIRST STEP. THE FIRST STEP * C* MUST BE DONE WITH THIS VALUE OF JSTART * C* SO THAT THE SUBROUTINE CAN INITIALIZE * C* ITSELF. * C* *1 TAKE A NEW STEP CONTINUING FROM THE LAST. * C* JSTART 16 SET TO NQ, THE CURRENT ORDER OF THE METHOD * C+ AT EXIT. NQ IS ALSO THE ORDER OF THE MAXIMUM * C* DERIVATIVE AVAILABLE. * C * PW A BLOCK OF AT LEAST N**2 FLOATING POINT LOCATIONS. * i ********* + **♦*** + ****♦*** + ♦ + **♦**** + * + + **** + ** + * + + + **♦ + + + + **** + + + .♦.** + ♦* DIMENSION V,VMAX(N>,SAVE,DY# 1 CSAVECN, 2) , IP(N) IRET - 1 KFLAO ■ 1 IF k JSTART. LE.0> GO TO 140 3U DO 110 I » 1 ,N DO lie J - 1 ,K SAVEC I , J) - V( I , J) HOLD - HNEW IF '. H . E Q . H L D ) GO TO 130 120 I RET 1 - 1 GO TO 7 5 I 3 6 N OLD = N Q TOLD = T IF (JSTART. 6T. 8) GO TO 250 GO TO 170 140 IF (JSTART.EQ.-I) GO TO 160 r* + ******* + + * + ** + * + *** + * + *** + * + * + *****# + ********** + + **** + ** + + + * + + + ** + + ■*•+ O ON THE FIR3T CALL, THE ORDER IS SET TO 1 AND THE INITIAL + C* DERIVATIVES PRE CALCULATED. * BR - 1.0 HO. - 1 C CALL DIFFUNCT, V, DV) CALL LINK< ' DIFFUN ' > DO 150 I - 1,N 1 50 V< I , 2) - DV< I )*H HNEW - H K - 2 GO TO 100 BO TO 100 C* REPEAT LAST STEP BV RESTORING SAVED INFORMATION. r »** + ***** + * *. * * **** + + ** + #*** + *** + *>c*** + **** + * + * + + ****** + * + **i('** + + 160 IF + * + * + + * + + * + ** + + + + + * + + + *** + * + + + *'C + + +* * UP TO 2 CORRECTOR ITERATIONS ARE TAKEN. CONVERGENCE IS TESTED BV + * REQUIRING THE L2 NORM OF CHANGES TO BE LESS THAN bND WHICH IS + + DEPENDENTONTHEERRORTESTCONSTANT. * THE SUM OF THE CORRECTIONS IS ACCUMULATED IN OR(I>. IT IS EQUAL TO THE K-TH DERIVATIVE H**K' i / , Hni' i» inLMLrunE PROPORTIONAL THE LOWEST POWER OF H PRESENT. => DVi. + + * + * + + C+ IF THERE HAS BEEN A CHANGE OF ORDER OR THERE HAS BEEN TROUBLE * C* WITH CONVERGENCE, PW IS RE-EVALUATED PRIOR TO STARTING THE * C * CORRECTOR ITERATION IN THE CASE OF STIFF METHODS. IWEVAL IS * C* THEN SET TO -1 AS AN INDICATOR THAT IT HAS BEEN DONE. * 310 DO 340 J ■= 1 , N F = V< J, 1 ) 330 340 R ■= EPS+AMAX1 *D 35 C + ****>!. + *****#*#******#*♦*••*****••*****************•** ********* 41*41 4141414141 C* ADD THE H>ENTITV MATRIX TO THE JACOBIAN AND DECOMPOSE INTO LU - PH * C#* + + + +**#****.*» + ** + ***#«i**4i*»*«*i*** + *****«i4l***** + fr* + ******* + **** + + * + ** + 296 DO 386 I ■ 1 H 380 PW< I*(N+1 )-N> - 1.0 ♦ PN-M> IUEVRL - -1 C 358 CALL SUB7( N, V, SAVE, C8AVE , VMAX , ERROR , PH, IP > 330 COLL L1NKC ' SUB7 ' > GO TO (448,410), JLNK C 410 CALL SU84(N, V. SAVE, CSAVE, VMAX, ERROR, PH, IP) 410 CALL LINK( 'SUB4' > GO TO <430, 448), JLNK 430 CONTINUE C 440 CPLL SUB3(N, V , SAVE , Ct AVE , VNAX, ERROR, PH, IP) 440 CPLL LINKC 'SUBS' ) GO TO )> JSTART - NQ 715 RETURN 720 IF(NQ.GT. 1 ) 80 TO 726 IF( ABS(H).LE.2.E8*HHIN> 80 TO 780 GO TO 550 725 R - H'HOLD DO 738 I - 1,N V( I , 1 > - SPVE( 1,1) 738 V(I,2) ■ SAVE 1 GO TO 170 740 KFLAG ■ -1 HNEH - H JSTART - NQ RETURN C 750 CALL SUB8 DIMENSION V(2,7>, DV(2) DV( 1 > - V(2, 1 ) DV<2> - -V( 1,1) RETURN END 36 Regular Overlays; INTRO vfth thq appropriate aubroutine name lotrodpcaa eaqft overly [NTHO .TITLE IN.r .GL06L N, RETURN .CSECT ODDR ODDR: .-. . CSECT BEGIN: ,SUB1 LIST! RET: J .END BEGIN I 15, SUB1 BR RET . word n, 0, 0,0,0, 0, e,e JSR 15, RETURN M c * * * + C* SE C* TV C* TH C* BE C* ST C* SC C**** c * * * * C* TH O c* c* c* c* c* c* c* OC sue COM 1 I 2ENQ DIM 1 ***,* T TH PE. IS S CPUS EP I OLIN * * * * 0(2 IF GO * * * * E FO CUfift ROUTINE HON 'CO DOUB, IN 2.ENQ3, ENSION ER ******* E COEFF CHECK ECTION E OF TH F IT HP G BEFOR ******* > * -1 . (NQ ,GT TO (221 ******* LLONING CV PERM SU61 (N, ns r<7>, EVPL, IRE BND.BR.D V,S ROR , P * * * * * * ** ICIE FOR SET E OR S NO E EX * * * * 8 . 6) ,222 COEFFIC ITTED BY NTS EXCE IWE DER T YE IT I *** * NQ- ,223 V,BAVE, JLNK, JS T, IRET1 EL, DELI 0VE(N,7 W(N,N), ******* THAT DE SSIVE YOL .8T CHANGE, T BEEN F IT HO ******* C60VE, TORT.K ,H,HOL ,D >, esov I P < N > ****** TERMIN RDER. .0 IF OND T DONE < S BEEN ****** YMOX, ERROR, PH, IP) ,KFLOa,NQ,NQOLD,NEWQ, EPS, T, TOLD, D,HNEW,HMfiX,HMIN, E,EUP,EDUN, ENQ1 , E OR SKIP TO FINAL ♦ COMPLETED - - 1 GO TO 230 ft< 1 > P<3) GO TO PC 1 > ft ( 3 ) P(4> GO TO 0< 1 ) ft ( 3 ) .0(4) ft < 5 > GO TO ft< 1 > •1/5,-1/50 *** * * * * * * 3/274J-83/274, -15/274, -1x274 /63,-15/36,-23/25 2,-3/232,-lx ***************************** 1 764 *********+******+********. -0, -0, 238 238 0(3) ft < 4 ) 0(5) ft<6> GO TO ft <: i > - fi<3> - ft < 4 > - 0(5) = 0(6) » ft<7> • RETURN END ft( -8 230 -0 -0 -0 -0 230 -0 -0 -8 -8 -0 230 -8 -0 -8 -8 -0 -8 809808888 6666666666666667 3333333333333333 8. 5454543454545453 1 > . 89898909098909891 . 480088088 . 708008800 .280 8 00000 . 020008080 . 437956204379562 . 821 1 6790321 1 6788 .3 102189781821898 . 05474452S54744526 . 8836496350364963384 .4861632653861 225 . 9206349206349286 . 4 1 * 666 66 6 666 66 6 7 . 0992863492863492 .81 1 984761 984761 9 .008566893424036282 37 SUBROUTINE 8UB2 < N , V , 8AVE , C8 AVE , VMAX , ERROR , PW , IP) COMMON •COM/' A, JUNK, JST ART , K. KFL AQ , NQ , NQOLD , NEWQ , EPS, T, TOLD, 1 IDOUB, IWEVAL, I RET, IRET1 , H, HOLD, HNEW, HMAX, HM I N, E , EUP, EDWN, ENQ1 , 2ENQ2,ENQ3,BND,BR,DEL,DEL1 , D DIMENSION V,CSAVE, VMAXCN), 1 ERROR,PM, IP DIMENSION PERTST<7,3> C* + + + + * + + + * + + + + + * + + + + + **m + ** + *<* + + + +-*** + + + + + + m + + + * + + * + **. + ± + + + ** + + 4,* + + + t, + 4, C* THE COEFFICIENTS IN PERT8T ARE USED IN 8ELECTINQ THE STEP AND * C* ORDER, THEREFORE ONLV ABOUT ONE PERCENT PCCURRCV IS NEEDED. * DATA PERTST/2.,4.5,7.333,10.42,13,7,17.15,1., 1 3. ,6. , 9. 1 67, 12.5, 15.98, 1 . , 1 . , 2 I ., 1 .,. 5, . 1 667, . 94U7, . 808333-' JLNK » 1 230 K - NCi+l IDOUB - K ENQ2 ■ . 5/FL0AT ENG3 - .5/FL0AT ENQ1 - 0.5/FLOAT EUP - *CPS>**2 E - *EPS>**2 EDWN - COMMON 'COM' A<7) , JLNK, JSTRRT, K, KFLAQ, NQ, NQOLD, NEWQ, EPS, T, TOLD, 1 IDOUB, IWEVAL, IRET.IRET1 , H , HOLD , HNEW , HMAX , HM I N , E , EUP , EDWN , ENQ 1 , 2ENQ2,ENQ3,BND,BR,DEL,DEL1 , D DIMENSION V(N,7),SAVE , IP C******+********. ***+******•*******************************♦.*** **♦+**••« C* THIS SECTION COMPUTES THE PREDICTED VALUES BV EFFECTIVELY * C* MULTIPLVINQ THE SAVED INFORMATION BV THE PASCAL TRIANGLE * C* MATRIX. * £••*******+*********•**+******•*******+*****+*******+*******+♦+*+*++*•** 250 T - T ♦ H DO 260 J - 2,K DO 268 J1 - J,K J2 ■ K - J1 ♦ J - 1 DO 260 I - 1 , N 260 V /COM/ A< 7 > JLNK, J8TART,K,KFLQG, NO, NQOLD, NEWQ, EPS, T, TOLD, , IWEVAL, I RET, 1 RET 1 , H , HOLD , HNEW , HMAX , HM I N , E , EUP , EDWN , ENQ 1 , 03, 0ND.BR, DEL, DELI , « u , M,VMAX .2) 0R< I EL TINUI . 2> DEL 1 .LE. BND) JLNK - 2 - V< I, 1 ) ♦ A<1 >*CSAVE - V< I , 2) - CSAVE< 1,1) ) - ERR0R< I ) ♦ CSAVE< 1,1) < C S A V E < I , 1 >/VMAXU)>**2 E JE BR - ANAX1 < . 9*BR, DEL/DEL1) 38 c * * + * C+ TH C+ PO O. TH C* MA C * IS C* ST C * * + + 440 468 470 480 SUBROUTI COMMON / t IDOUB, 2ENQ2, ENQ DIMENS 10 1 ********* E CORRECT SSIBILITI IS IS EIT TRIX PW H TAKEN. EP IS RED ********* GO TO (4 T = TOLD I T < QBS C H IF< IHEVP IWEVPL - I RET 1 - JLNK - 1 GO TO 33 KFLPG ■ NO - NGO DO 480 I DO 460 V.LE L. NE 2 2 6 -3 Lt> - 1 J - J) - NO UB5,C8flVECN,2),VMPK, PN, IP ************»**********•***********+*»**************** TERPTJON FAILED TO CONVERSE IN 2 TRIES. VARIOUS * RE CHECKED FOR. IF H 19 ALREADV HMIN AND * ADAMS METHOD OR THE STIFF METHOD IN WHICH THE * LREADV BEEN RE -EVALUATED , A NO CONVERGENCE EXIT * WISE THE MATRIX PW IS RE-EVALUATED AND'OR THE * TO TRV AND GET CONVERGENCE. * a***************************************************** 90,470, 470), JLNK .HMIN .AND. < IWEVAL-1 ) .LT. -1 > GO TO 463 .0) H - H*0,2SE0 ,N 1 ,K SOVC C* THE CORRECTOR CONVERSED AND CONTROL 16 PASSED TO 9TATEMENT 320 • C* IF THE ERROR TE6T 19 O.K., AND TO 540 OTHERWISE. * C* IF THE STEP IS O.K. IT 19 ACCEPTED. IF IDOUB HAS BEEN REDUCED * C* TO ONE, A TEST IS MADE TO 8EE IF TME 9TEP CAN BE INCREASED * C* PT THE CURRENT ORDER OR BV GOING TO ONE HIGHER OR ONE LOWER. • C* SUCH H CHAN8E IS ONLV MADE IF THE STEP CAN BE INCREASED BV AT * C* LEAST 1.1. IF NO CHANGE IB POSSIBLE IDOUB IS SET TO 8 TO ♦ C* PREVENT FUTHER TESTING FOR 9 STEPS. * C* IF P CHAN6E IS POSSIBLE, IT IS MADE AND IDOUB IS SET TO • C* NQ ♦ 1 TO PREVENT FURTHER TESING FOR THAT NUMBER OF STEPS. • C* IF THE ERROR HAS TOO LARGE, THE OPTIMUM STEP SIZE FOR THIS OR * C* LOHER ORDER IS COMPUTED, AND THE 9TEP RETRIED. IF IT SHOULD * C* FAIL TWICE MORE IT IS AN INDICATION THAT THE DERIVOTIVES THAT * C* HAVE ACCUMULATED IN THE V ARRAV HAVE ERRORS OF THE WRONG ORDER * C* SO THE FIRST DERIVATIVES ARE RECOMPUTED AND THE ORDER IS SET * C* TO 1 . * 1 ,N (. ERROR* I >/VMAX< I > >**2 496 D » 0.0 DO 580 I 508 D ■ D IWEVPL - IF(D.GT.E) JLNK * 5 IF< JLNK. EQ. 5) GO TO 536 IF JLNK = 4 RETURN END j> ♦ r< j:>*error< i > ) JLNK-3 GO TO 536 535 39 SUBROUTINE 6U86 COMMON 'COM' fl<7>.JLNK. J«T«RT , K , KFLhO, NQ , NQOLD , NEWQ , EPS, T, TOLD, 1 IDOUB, IWEVfiL, I RET, IRET1 , H, HOLD , HNEW , HMfiX , HM I N , E,EUP,EDWN,ENQ1, 2ENQ2,ENU3,BND,BR,DEL, DELI , D , *■ u n u u l v dc u i v i u c v hi vn C+ ONE HIGHER RESPECT I VELV . GO TO (680, 358, 640) , JLNK 550 PR2 - >80 TO 576 r\ ^ a a 560 = 0.0 DO 560 I - 1 ,N ♦ <>'VMfiX>**2 PR3 « >«"»2 PR1 - 00 TO 650 IF 80 TO 660 R - 1 .©^flMfiXI NEWQ * NQ - 1 610 IDOUB - 8 IF ( ) > 60 TO 694 IK - ERROR* I >*A GO TO 695 NO ■ NEMQ JLNK » 3 GO TO 695 IF 00 TO 600 NEWQ - NO R - 1 . 0'HMfiX1 GO TO 6 10 R - 1.8'HMAX1tPR3,1.E-4) NEWQ - NQ ♦ 1 '3 TO 618 IRET - 2 R ■ AMIN1 H - H*R KNEW - H IF «:Nw. EO.NEWO) GO TO 688 NQ - NEWQ JLNK - 3 GO TO 695 R1 » 1.0 DO * 9 8 J - 2 , K R1 - R1 *R DO 6 9 8 I - 1 , N V < I , J > - V < I , J > * R 1 IDOUB = K JLNK - 1 RETURN END ++**+++ OUTPUT. fl< I , J P < I , J IP(K> IP USE 'SO DETERM< IF IPCN iRDER OF MPTRIX DECLARED DIMENSION OF PRftftV P. IPTRIX TO BE TRIANGULRRJZED. I LE.J » UPPER TRIANGULAR FACTOR, U. I.QT.J - MULTIPLIERS • LOWER TRIPNQULPR FACTOR, I-L, K.LT.N » INPEX OF K-TH PIVOT ROW, <-1 >** - IP*fl<1,1>*fl<2,2)*.,,*fl, > «= 0, A IS SINGULAR, SOLVE WILL DIVIDE BV ZERO, f 01 TO 5 ■IP 350 360 DIMENSION A ,B< 1 > , IP, V< 1 > GO TO ■; 1 1 , 3 5 ) , J L N K J L N K = 2 N D I M = N I P < N > - 1 DO 6 K=1 , N I F < K . E 6 . N > GO KP1 - K+1 U m U DO 1 I-KP1,N IF.QT.PBS>> CONTINUE IP - M IF(M.NE.K) IP T - fl P < M , K ) - P t K , K > fl < K , K > •» T I F < T . E Q . > GO TO DO 2 I = K P 1 P< I , K> ■ DO 4 J-KP1 , N T - fl < M , J > fl ♦ R(| ,K>*T CONTINUE IF.EQ.0 CONTINUE IF< IP . EO. 0> JLNK ■ 1 IF< JLNK. EQ. 1 > 360 Ixl.N M-l N -P< I , K> / T ' GO TO 4 P< I, J> IP(N> IF IFC JLNK. EO. 1 ) GO TO 10 DO B< I > JLNK V< I+N> - BC I >*H SOLUTION OF LINEAR SVSTEM, A*X - B. INPUT. . . N - ORDER OF MPTRIX. NDIM - DECLARED DIMENSION OF PRRPV P. fl =• TRIANGULARIZED MATRIX OBTAINED FROM 6 - RIGHT HAND SIDE VECTOR. IP = PIVOT VECTOR OBTAINED FROM 'DECOMP OUTPUT. . . B - SOLUTION VECTOR, X. I F < N . E 6 . 1 > GO TO « NM1 = N-1 DO 7 K - 1 , N M 1 KP1 = K + 1 M - IP T = B < M > K)*T K = K M 1 + 1 B < K > = B < K > / fl ( K , K > T = - 6 < K > 00 8 1=1 ,KM1 * B < I > - B < I > ♦ A < I , K "> B < 1 > = B < 1 ) / P < 1 . 1 > 10 RETURN END DECOMP B(M ■■> - B B < K ) = T DO 7 I-KP1 , N B< I •> = BCD + P< I DO '8 KB-1 , NM1 KM1 = : N - KB i *T ill C* TH C * TO C + *• * * 759 76S 770 76* 783 SUBROUTINE 6UB6 < N , V , SAVE , CS AVE , VM AH , ERROR , PW , IP > COMMON -'COM/ A < 7 > , JLNK , J8TART , K , KFLfiS , NO. , NQOLD , NEWO , EPS , T , TOLD , 1 IC)OUB,IWEVPL I IRET > !RETl,M,H01.D,HNEW,HMPX ) HMIN,E I EUP ) Et'WN,ENi?1, 2ENQ2,ENQ3,BND,BR, DEL, DEL 1 , t> DIMENSION V(N,7),S«VE, IP *V++++>**+4V + + + * * + *** + * + **-K* + ****** + *** + **t.+ + + + + **** + 4'* + + + + + + + + +-l<* + IS SECTION SCALES ALL VARIABLES CONNECTED WITH H AND RETURNS * THE ENTERING SECTION. * H = RMfiXI R1 = 1.0 DO 7 - SAVE< I , J>*R1 DO 7 7 I - 1 , N V < I , 1) » S R V E < I , 1 > IDOUB - K JLNK » I R E T 1 GO TO 7 8 5 K.KLAG - -4 JLNK - 4 RETURN END Form AEC-427 (6/68) AECM 3201 U.S. ATOMIC ENERGY COMMISSION UNIVERSITY-TYPE CONTRACTOR'S RECOMMENDATION FOR DISPOSITION OF SCIENTIFIC AND TECHNICAL DOCUMENT ( See Instructions on Reverse Side ) 1. AEC REPORT NO. COO-ll|69-0215 2. TITLE NUMERICAL SYSTEMS ON A MINICOMPUTER Roy Leonard Brown, Jr. 3. TYPE OF DOCUMENT (Check one): P^l a. Scientific and technical report | I b. Conference paper not to be published in a journal: Title of conference Date of conference Exact location of conference. Sponsoring organization □ c. Other (Specify) 4. RECOMMENDED ANNOUNCEMENT AND DISTRIBUTION (Check one): [x] a. AEC's normal announcement and distribution procedures may be followed. Q b. Make available only within AEC and to AEC contractors and other U.S. Government agencies and their contractors. ~] c. Make no announcement or distrubution. 5. REASON FOR RECOMMENDED RESTRICTIONS: 6. SUBMITTED BY: NAME AND POSITION (Please print or type) C. W. Gear Professor and Principal Investigator Organization Department of Computer Science University of Illinois Urbana, Illinois 6l801 Signature ^J^hMl^ ^a^-. Date February 1973 FOR AEC USE ONLY AEC CONTRACT ADMINISTRATOR'S COMMENTS, IF ANY, ON ABOVE ANNOUNCEMENT AND DISTRIBUTION RECOMMENDATION: PATENT CLEARANCE: [_) a. AEC patent clearance has been granted by responsible AEC patent group. LJ b. Report has been sent to responsible AEC patent group for clearance. LI c. Patent clearance not required. IBLIOGRAPHIC DATA 1EET 1. Report No. UIUCDCS-R-73-555 3. Recipient's Accession No. Title and Subtitle NUMERICAL SYSTEMS ON A MINICOMPUTER 5- Report Date February 1973 6. Author(s) R oy Leonard Brovn , Jr . 8. Performing Organization Rept. No - COO-1U69-0215 Performing Organization Name and Address Department of Computer Science University of Illinois Urban a, Illinois 6l801 10. Project/Task/Work Unit No. US AEC AT (11-1)11+69 11. Contract /Grant No. Sponsoring Organization Name and Address US AEC Chicago Operations Office 9800 South Cass Avenue Argonne, Illinois 60U39 13. Type of Report & Period Covered Thesis Research 14. Supplementary Notes Abstracts This thesis defines the concept of a numerical system for a minicomputer and provides a description of the software and computer system configuration necessary to implement such a system. A procedure for creating a numerical system from a FORTRAN program is developed and an example is presented. The reader should have some knowledge of FORTRAN and minicomputer operating systems, PAL assembly language for PDP-11 1 Key Words and Document Analysis. 17a. Descriptors minicomputer system numerical system overlay buffer overlay partition of FORTRAN program overlay facility 17l Identifiers/Open-Ended Terms 9 COSATI Field/Group Plvailability Statement unlimited distribution 19. Security Class (This Report) UNCLASSIFIED 20. Security Class (This Page UNCLASSIFIED 21. No. of Pages kQ 22. Price P °' NTIS-35 (10-70) USCOMM-DC 40329-P71 OCT 2 4 19/3