home *** CD-ROM | disk | FTP | other *** search
/ NetNews Usenet Archive 1992 #18 / NN_1992_18.iso / spool / comp / sys / intel / 1535 < prev    next >
Encoding:
Internet Message Format  |  1992-08-19  |  42.7 KB
  1. Path: sparky!uunet!kithrup!hoptoad!pacbell.com!mips!swrinde!elroy.jpl.nasa.gov!usc!sol.ctr.columbia.edu!ira.uka.de!uka!uka!news
  2. From: S_JUFFA@iravcl.ira.uka.de (|S| Norbert Juffa)
  3. Newsgroups: comp.sys.intel
  4. Subject: What you always wanted to know about math coprocessors for 80x86 4/4
  5. Message-ID: <16tnskINNcs1@iraul1.ira.uka.de>
  6. Date: 19 Aug 92 15:03:48 GMT
  7. Organization: University of Karlsruhe (FRG) - Informatik Rechnerabt.
  8. Lines: 903
  9. NNTP-Posting-Host: irav1.ira.uka.de
  10. X-News-Reader: VMS NEWS 1.23
  11.  
  12.  
  13.          References
  14.  
  15.          [1]  Schnurer, G.: Zahlenknacker im Vormarsch.
  16.               c't 1992, Heft 4, Seiten 170-186
  17.          [2]  Curnow, H.J.; Wichmann, B.A.: A synthetic benchmark.
  18.               Computer Journal, Vol. 19, No. 1, 1976, pp. 43-49
  19.          [3]  Wichmann, B.A.: Validation code for the Whetstone benchmark.
  20.               NPL Report DITC 107/88, National Physics Laboratory, UK,
  21.               March 1988
  22.          [4]  Curnow, H.J.: Wither Whetstone? The Synthetic Benchmark after
  23.               15 Years.
  24.               In: Aad van der Steen (ed.): Evaluating Supercomputers.
  25.               London: Chapman and Hall 1990
  26.          [5]  Dongarra, J.J.: The Linpack Benchmark: An Explanation.
  27.               In: Aad van der Steen (ed.): Evaluating Supercomputers.
  28.               London: Chapman and Hall 1990
  29.          [6]  Dongarra, J.J.: Performance of Various Computers Using Standard
  30.               Linear Equations Software.
  31.               Report CS-89-85, Computer Science Department, University of
  32.               Tennessee, March 11, 1992
  33.          [7]  Huth, N.: Dichtung und Wahrheit oder Datenblatt und Test.
  34.               Design & Elektronik 1990, Heft 13, Seiten 105-110
  35.          [8]  Ungerer, B.: Sockelfolger.
  36.               c't 1990, Heft 4, Seiten 162-163
  37.          [9]  Coonen, J.T.: Contributions to a Proposed Standard for Binary
  38.               Floating-Point Arithmetic
  39.               Ph.D. thesis, University of California, Berkeley, 1984
  40.          [10] IEEE: IEEE Standard for Binary Floating-Point Arithmetic.
  41.               SIGPLAN Notices, Vol. 22, No. 2, 1985, pp. 9-25
  42.          [11] IEEE Standard for Binary Floating-Point Arithmetic.
  43.               ANSI/IEEE Std 754-1985.
  44.               New York, NY: Institute of Electrical and Electronics
  45.               Engineers 1985
  46.          [12] FasMath 83D87 Compatibility Report. Cyrix Corporation, Nov. 1989
  47.               Order No. B2004
  48.          [13] FasMath 83D87 Accuracy Report. Cyrix Corporation, July 1990
  49.               Order No. B2002
  50.          [14] FasMath 83D87 Benchmark Report. Cyrix Corporation, June 1990
  51.               Order No. B2004
  52.          [15] FasMath 83D87 User's Manual. Cyrix Corporation, June 1990
  53.               Order No. L2001-003
  54.          [16] Brent, R.P.: A FORTRAN multiple-precision arithmetic package.
  55.               ACM Transactions on Mathematical Software, Vol. 4, No. 1,
  56.               March 1978, pp. 57-70
  57.          [17] 387DX User's Manual, Programmer's Reference. Intel Corporation,
  58.               1989
  59.               Order No. 231917-002
  60.          [18] Volder, J.E.: The CORDIC Trigonometric Computing Technique.
  61.               IRE Transactions on Electronic Computers, Vol. EC-8, No. 5,
  62.               September 1959, pp. 330-334
  63.          [19] Walther, J.S.: A unified algorithm for elementary functions.
  64.               AFIPS Conference Proceedings, Vol. 38, SJCC 1971, pp. 379-385
  65.          [20] Esser, R.; Kremer, F.; Schmidt, W.G.: Testrechnungen auf der
  66.               IBM 3090E mit Vektoreinrichtung.
  67.               Arbeitsbericht RRZK-8803, Regionales Rechenzentrum an der
  68.               Universit"at zu Kln, Februar 1988
  69.          [21] McMahon, H.H.: The Livermore Fortran Kernels: A test of the
  70.               numerical performance range.
  71.               Technical Report UCRL-53745, Lawrence Livermore National
  72.               Laboratory, USA, December 1986
  73.          [22] Nave, R.: Implementation of Transcendental Functions on a Numerics
  74.               Processor.
  75.               Microprocessing and Microprogramming, Vol. 11, No. 3-4,
  76.               March-April 1983, pp. 221-225
  77.          [23] Yuen, A.K.: Intel's Floating-Point Processors.
  78.               Electro/88 Conference Record, Boston, MA, USA, 10-12 May 1988,
  79.               pp. 48/5-1 - 48/5-7
  80.          [24] Stiller, A.; Ungerer, B.: Ausgerechnet.
  81.               c't 1990, Heft 1, Seiten 90-92
  82.          [25] Rosch, W.L.: Handfeste Hilfe oder Seifenblase?
  83.               PC Professionell, Juni 1991, Seiten 214-237
  84.          [26] Intel 80286 Hardware Reference Manual. Intel Corporation, 1987
  85.               Order No.210760-002
  86.          [27] AMD 80C287 80-bit CMOS Numeric Processor. Advanced Micro Devices,
  87.               June 1989
  88.               Order No. 11671B/0
  89.          [28] Intel RapidCAD(tm) Engineering CoProcessor Performance Brief.
  90.               Intel Corporation, 1992
  91.          [29] i486(tm) Microprocessor Performance Report. Intel Corporation,
  92.               April 1990
  93.               Order No. 240734-001
  94.          [30] Intel486(tm) DX2 Microprocessor Performance Brief. Intel
  95.               Corporation, March 1992
  96.               Order No. 241254-001
  97.          [31] Abacus 3167 Floating-Point Coprocessor Data Book. Weitek
  98.               Corporation, July 1990
  99.               DOC No. 9030
  100.          [32] WTL 4167 Floating-Point Coprocessor Data Book. Weitek
  101.               Corporation, July 1989
  102.               DOC No. 8943
  103.          [33] Abacus Software Designer's Guide. Weitek Corporation,
  104.               September 1989
  105.               DOC No. 8967
  106.          [34] Stiller, A.: Cache & Carry.
  107.               c't 1992, Heft 6, Seiten 118-130
  108.          [35] Stiller, A.: Cache & Carry, Teil 2.
  109.               c't 1992, Heft 7, Seiten 28-34
  110.          [36] Palmer, J.F.; Morse, S.P.: Die mathematischen Grundlagen der
  111.               Numerik-Prozessoren 8087/80287.
  112.               Mnchen: tewi 1985
  113.          [37] 80C187 80-bit Math Coprocessor Data Sheet. Intel Corporation,
  114.               September 1989
  115.               Order No. 270640-003
  116.          [38] IIT-2C87 80-bit Numeric Co-Processor Data Sheet. IIT, May 1990
  117.          [39] Engineering note 4x4 matrix multiply transformation. IIT, 1989
  118.          [40] Tscheuschner, E.: 4 mal 4 auf einen Streich.
  119.               c't 1990, Heft 3, Seiten 266-276
  120.          [41] Goldberg, D.: Computer Arithmetic.
  121.               In: Hennessy, J.L.; Patterson, D.A.: Computer Architecture A
  122.               Quantitative Approach. San Mateo, CA: Morgan Kaufmann 1990
  123.          [42] 8087 Math Coprocessor Data Sheet. Intel Corporation, October 1989,
  124.               Order No. 205835-007
  125.          [43] 8086/8088 User's Manual, Programmer's and Hardware Reference.
  126.               Intel Corporation, 1989
  127.               Order No. 240487-001
  128.          [44] 80286 and 80287 Programmer's Reference Manual. Intel Corporation,
  129.               1987
  130.               Order No. 210498-005
  131.          [45] 80287XL/XLT CHMOS III Math Coprocessor Data Sheet. Intel
  132.               Corporation, May 1990
  133.               Order No. 290376-001
  134.          [46] Cyrix FasMath(tm) 82S87 Coprocessor Data Sheet. Cyrix Coporation,
  135.               1991
  136.               Document 94018-00 Rev. 1.0
  137.          [47] IIT-3C87 80-bit Numeric Co-Processor Data Sheet. IIT, May 1990
  138.          [48] 486(tm)SX(tm) Microprocessor/ 487(tm)SX(tm) Math CoProcessor
  139.               Data Sheet. Intel Corporation, April 1991.
  140.               Order No. 240950-001
  141.          [49] Schnurer, G.: Die gro"se Verlade.
  142.               c't 1991, Heft 7, Seiten 55-57
  143.          [50] Schnurer, G.: Eine 4 f"ur alle.
  144.               c't 1991, Heft 6, Seite 25
  145.          [51] Intel486(tm)DX Microprocessor Data Book. Intel Corporation,
  146.               June 1991
  147.               Order No. 240440-004
  148.          [52] i486(tm) Microprocessor Hardware Reference Manual. Intel
  149.               Corporation, 1990
  150.               Order No. 240552-001
  151.          [53] i486(tm) Microprocessor Programmer's Reference Manual. Intel
  152.               Corporation, 1990
  153.               Order No. 240486-001
  154.          [54] Ungerer, B.: Kalte H"ute.
  155.               c't 1992, Heft 8, Seiten 140-144
  156.          [55] Ungerer, B.: Hei"se Sache.
  157.               c't 1991, Heft 4, Seiten 104-108
  158.          [56] Rosch, W.L.: Handfeste Hilfe oder Seifenblase?
  159.               PC Profesionell, Juni 1991, Seiten 214-237
  160.          [57] Niederkr"uger, W.: Lebendige Vergangenheit.
  161.               c't 1990, Heft 12, Seiten 114-116
  162.          [58] ULSI Math*Co Advanced Math Coprocessor Technical Specification.
  163.               ULSI System, 5/92, Rev. E
  164.          [59] 387(tm)DX Math CoProcessor Data Sheet. Intel Corporation,
  165.               September 1990.
  166.               Order No. 240448-003
  167.          [60] 387(tm) Numerics Coprocessor Extension Data Sheet. Intel
  168.               Corporation, February 1989.
  169.               Order No. 231920-005
  170.          [61] Koren, I.; Zinaty, O.: Evaluating Elementary Functions in a
  171.               Numerical Coprocessor Based on Rational Approximations.
  172.               IEEE Transactions on Computers, Vol. C-39, No. 8, August 1990,
  173.               pp. 1030-1037
  174.          [62] 387(tm) SX Math CoProcessor Data Sheet. Intel Corporation,
  175.               November 1989
  176.               Order No. 240225-005
  177.          [63] Frenkel, G.: Coprocessors Speed Numeric Operations.
  178.               PC-Week, August 27, 1990
  179.          [64] Schnurer, G.; Stiller, A.: Auto-Matt.
  180.               c't 1991, Heft 10, Seiten 94-96
  181.          [65] Grehan, R.: FPU Face-Off.
  182.               Byte, November 1990, pp. 194-200
  183.          [66] Tang, P.T.P.: Testing Computer Arithmetic by Elementary Number
  184.               Theory. Preprint MCS-P84-0889, Mathematics and Computer Science
  185.               Division, Argonne National Laboratory, August 1989
  186.          [67] Ferguson, W.E.: Selecting math coprocessors.
  187.               IEEE Spectrum, July 1991, pp. 38-41
  188.          [68] Schnabel, J.: Viermal 387.
  189.               Computer Pers"onlich 1991, Heft 22, Seiten 153-156
  190.          [69] Hofmann, J.: Starke Rechenknechte.
  191.               mc 1990, Heft 7, Seiten 64-67
  192.          [70] Woerrlein, H.; Hinnenberg, R.: Die Lust an der Power.
  193.               Computer Live 1991, Heft 10, Seiten 138-149
  194.  
  195.  
  196.  
  197.          Manufacturer's addresses
  198.  
  199.          Intel Corporation
  200.          3065 Bowers Avenue
  201.          Santa Clara, CA 95051
  202.          USA
  203.  
  204.          IIT Integrated Information Technology, Inc.
  205.          2540 Mission College Blvd.
  206.          Santa Clara, CA 95054
  207.          USA
  208.  
  209.          ULSI Systems, Inc.
  210.          58 Daggett Drive
  211.          San Jose, CA 95134
  212.          USA
  213.  
  214.          Chips & Technologies, Inc.
  215.          3050 Zanker Road
  216.          San Jose, CA 95134
  217.          USA
  218.  
  219.          Weitek Corporation
  220.          1060 East Arques Avenue
  221.          Sunnyvale, CA 94086
  222.          USA
  223.  
  224.          AMD Advanced Microdevices, Inc.
  225.          901 Thompson Place
  226.          P.O.B. 3453
  227.          Sunnyvale, CA 94088-3453
  228.          USA
  229.  
  230.          Cyrix Corporation
  231.          P.O.B. 850118
  232.          Richardson, TX 75085
  233.          USA
  234.  
  235.  
  236.          Appendix A
  237.  
  238.  
  239.          {$N+,E+}
  240.          PROGRAM PCtrl;
  241.  
  242.          VAR B,c: EXTENDED;
  243.              Precision, L: WORD;
  244.  
  245.          PROCEDURE SetPrecisionControl (Precision: WORD);
  246.          (* This procedure sets the internal precision of the NDP. Available *)
  247.          (* precision values:  0  -  24 bits (SINGLE)                        *)
  248.          (*                    1  -  n.a. (mapped to single)                 *)
  249.          (*                    2  -  53 bits (DOUBLE)                        *)
  250.          (*                    3  -  64 bits (EXTENDED)                      *)
  251.  
  252.          VAR CtrlWord: WORD;
  253.  
  254.          BEGIN {SetPrecisionCtrl}
  255.             IF Precision = 1 THEN
  256.                Precision := 0;
  257.             Precision := Precision SHL 8; { make mask for PC field in ctrl word}
  258.             ASM
  259.                FSTCW    [CtrlWord]        { store NDP control word }
  260.                MOV      AX, [CtrlWord]    { load control word into CPU }
  261.                AND      AX, 0FCFFh        { mask out precision control field }
  262.                OR       AX, [Precision]   { set desired precision in PC field }
  263.                MOV      [CtrlWord], AX    { store new control word }
  264.                FLDCW    [CtrlWord]        { set new precision control in NDP }
  265.             END;
  266.          END; {SetPrecisionCtrl}
  267.  
  268.          BEGIN {main}
  269.             FOR Precision := 1 TO 3 DO BEGIN
  270.                B := 1.2345678901234567890;
  271.                SetPrecisionControl (Precision);
  272.                FOR L := 1 TO 20 DO BEGIN
  273.                   B := Sqrt (B);
  274.                END;
  275.                FOR L := 1 TO 20 DO BEGIN
  276.                   B := B*B;
  277.                END;
  278.                SetPrecisionControl (3);   { full precision for printout }
  279.                WriteLn (Precision, B:28);
  280.             END;
  281.          END.
  282.  
  283.  
  284.          +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  285.  
  286.          {$N+,E+}
  287.          PROGRAM RCtrl;
  288.  
  289.          VAR B,c: EXTENDED;
  290.              RoundingMode, L: WORD;
  291.  
  292.  
  293.          PROCEDURE SetRoundingMode (RCMode: WORD);
  294.          (* This procedure selects one of four available rounding modes *)
  295.          (* 0  -  Round to nearest (default)                            *)
  296.          (* 1  -  Round down (towards negative infinity)                *)
  297.          (* 2  -  Round up (towards positive infinity)                  *)
  298.          (* 3  -  Chop (truncate, round towards zero)                   *)
  299.  
  300.          VAR CtrlWord: WORD;
  301.  
  302.          BEGIN
  303.             RCMode := RCMode SHL 10;      { make mask for RC field in control word}
  304.             ASM
  305.                FSTCW    [CtrlWord]        { store NDP control word }
  306.                MOV      AX, [CtrlWord]    { load control word into CPU }
  307.                AND      AX, 0F3FFh        { mask out rounding control field }
  308.                OR       AX, [RCMode]      { set desired precision in RC field }
  309.                MOV      [CtrlWord], AX    { store new control word }
  310.                FLDCW    [CtrlWord]        { set new rounding control in NDP }
  311.             END;
  312.          END;
  313.  
  314.          BEGIN
  315.             FOR RoundingMode := 0 TO 3 DO BEGIN
  316.                B := 1.2345678901234567890e100;
  317.                SetRoundingMode (RoundingMode);
  318.                FOR L := 1 TO 51 DO BEGIN
  319.                   B := Sqrt (B);
  320.                END;
  321.                   FOR L := 1 TO 51 DO BEGIN
  322.                   B := -B*B;
  323.                END;
  324.                SetRoundingMode (0);        { round to nearest for printout }
  325.                WriteLn (RoundingMode, B:28);
  326.             END;
  327.          END.
  328.  
  329.  
  330.          +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  331.  
  332.          {$N+,E+}
  333.  
  334.          PROGRAM DenormTs;
  335.  
  336.          VAR E: EXTENDED;
  337.              D: DOUBLE;
  338.              S: SINGLE;
  339.  
  340.          BEGIN
  341.             WriteLn ('Testing support and printing of denormals');
  342.             WriteLn;
  343.             Write ('Coprocessor is: ');
  344.             CASE Test8087 OF
  345.                0: WriteLn ('Emulator');
  346.                1: WriteLn ('8087 or compatible');
  347.                2: WriteLn ('80287 or compatible');
  348.                3: WriteLn ('80387 or compatible');
  349.             END;
  350.             WriteLn;
  351.             S := 1.18e-38;
  352.             S := S * 3.90625e-3;
  353.             IF S = 0 THEN
  354.                WriteLn ('SINGLE denormals not supported')
  355.             ELSE BEGIN
  356.                WriteLn ('SINGLE denormals supported');
  357.                WriteLn ('SINGLE denormal prints as:   ', S);
  358.                WriteLn ('Denormal should be printed as 4.60943...E-0041');
  359.             END;
  360.             WriteLn;
  361.             D := 2.24e-308;
  362.             D := D * 3.90625e-3;
  363.             IF D = 0 THEN
  364.                WriteLn ('DOUBLE denormals not supported')
  365.             ELSE BEGIN
  366.                WriteLn ('DOUBLE denormals supported');
  367.                WriteLn ('DOUBLE denormal prints as:   ', D);
  368.                WriteLn ('Denormal should be printed as 8.75...E-0311');
  369.             END;
  370.             WriteLn;
  371.             E := 3.37e-4932;
  372.             E := E * 3.90625e-3;
  373.             IF E = 0 THEN
  374.                WriteLn ('EXTENDED denormals not supported')
  375.             ELSE BEGIN
  376.                WriteLn ('EXTENDED denormals supported');
  377.                WriteLn ('EXTENDED denormal prints as: ', E);
  378.                WriteLn ('Denormal should be printed as 1.3164...E-4934');
  379.             END;
  380.          END.
  381.  
  382.  
  383.          Appendix B
  384.  
  385.  
  386.          ; FILE: APFELM4.ASM
  387.          ; assemble with MASM /e APFELM4 or TASM /e APFELM4
  388.  
  389.  
  390.          CODE        SEGMENT BYTE PUBLIC 'CODE'
  391.                      ASSUME  CS: CODE
  392.  
  393.                      PAGE    ,120
  394.  
  395.                      PUBLIC  APPLE87;
  396.  
  397.          APPLE87     PROC    NEAR
  398.                      PUSH    BP                  ; save caller's base pointer
  399.                      MOV     BP, SP              ; make new frame pointer
  400.                      PUSH    DS                  ; save caller's data segment
  401.                      PUSH    SI                  ; save register
  402.                      PUSH    DI                  ;  variables
  403.                      LDS     BX, [BP+04]         ; pointer to parameter record
  404.                      FINIT                       ; init 80x87          FSP->R0
  405.                      FILD   WORD  PTR [BX+02]    ; maxrad              FSP->R7
  406.                      FLD    QWORD PTR [BX+08]    ; qmax                FSP->R6
  407.                      FSUB   QWORD PTR [BX+16]    ; qmax-qmin           FSP->R6
  408.                      DEC    WORD  PTR [BX+04]    ; ymax-1
  409.                      FIDIV  WORD  PTR [BX+04]    ; (qmax-qmin)/(ymax-1)FSP->R6
  410.                      FSTP   QWORD PTR [BX+16]    ; save delta_q        FSP->R7
  411.                      FLD    QWORD PTR [BX+24]    ; pmax                FSP->R6
  412.                      FSUB   QWORD PTR [BX+32]    ; pmax-pmin           FSP->R6
  413.                      DEC    WORD  PTR [BX+06]    ; xmax-1
  414.                      FIDIV  WORD  PTR [BX+06]    ; delta_p             FSP->R6
  415.                      MOV    AX, [BX]             ; save maxiter,[BX] needed for
  416.                      MOV    [BX+2], AX           ;  80x87 status now
  417.                      XOR    BP, BP               ; y=0
  418.                      FLD    QWORD PTR [BX+08]    ; qmax                FSP->R5
  419.                      CMP    WORD  PTR [BX+40], 0 ; fast mode on 8087 desired ?
  420.                      JE     yloop                ; no, normal mode
  421.                      FSTCW  [BX]                 ; save NDP control word
  422.                      AND    WORD PTR [BX], 0FCFFh; set PCTRL = single precision
  423.                      FLDCW  [BX]                 ; get back NDP control word
  424.          yloop:      XOR    DI, DI               ; x=0
  425.                      FLD    QWORD PTR [BX+32]    ; pmin                FSP->R4
  426.          xloop:      FLDZ                        ; j**2= 0             FSP->R3
  427.                      FLDZ                        ; 2ij = 0             FSP->R2
  428.                      FLDZ                        ; i**2= 0             FSP->R1
  429.                      MOV    CX, [BX+2]           ; maxiter
  430.                      MOV    DL, 41h              ; mask for C0 and C3 cond.bits
  431.          iteration:  FSUB   ST, ST(2)            ; i**2-j**2           FSP->R1
  432.                      FADD   ST, ST(3)            ; i**2-j**2+p = i     FSP->R1
  433.                      FLD    ST(0)                ; duplicate i         FSP->R0
  434.                      FMUL   ST(1), ST            ; i**2                FSP->R0
  435.                      FADD   ST, ST(0)            ; 2i                  FSP->R0
  436.                      FXCH   ST(2)                ; 2*i*j               FSP->R0
  437.                      FADD   ST, ST(5)            ; 2*i*j+q = j         FSP->R0
  438.                      FMUL   ST(2), ST            ; 2*i*j               FSP->R0
  439.                      FMUL   ST, ST(0)            ; j**2                FSP->R0
  440.                      FST    ST(3)                ; save j**2           FSP->R0
  441.                      FADD   ST, ST(1)            ; i**2+j**2           FSP->R0
  442.                      FCOMP  ST(7)                ; i**2+j**2 > maxrad? FSP->R1
  443.                      FSTSW  [BX]                 ; save 80x87 cond.codeFSP->R1
  444.                      TEST   BYTE PTR [BX+1], DL  ; test carry and zero flags
  445.                      LOOPNZ iteration            ; until maxiter if not diverg.
  446.                      MOV    DX, CX               ; number of loops executed
  447.                      NEG    CX                   ; carry set if CX <> 0
  448.                      ADC    DX, 0                ; adjust DX if no. of loops<>0
  449.  
  450.                      ; plot point here (DI = X, BP = y, DX has the color)
  451.  
  452.                      FSTP   ST(0)                ; pop i**2            FSP->R2
  453.                      FSTP   ST(0)                ; pop 2ij             FSP->R3
  454.                      FSTP   ST(0)                ; pop j**2            FSP->R4
  455.                      FADD   ST,ST(2)             ; p=p+delta_p         FSP->R4
  456.                      INC    DI                   ; x:=x+1
  457.                      CMP    DI, [BX+6]           ; x > xmax ?
  458.                      JBE    xloop                ; no, continue on same line
  459.                      FSTP   ST(0)                ; pop p               FSP->R5
  460.                      FSUB   QWORD PTR [BX+16]    ; q=q-delta_q         FSP->R5
  461.                      INC    BP                   ; y:=y+1
  462.                      CMP    BP, [BX+4]           ; y > ymax ?
  463.                      JBE    yloop                ; no, picture not done yet
  464.  
  465.          groesser:   POP    DI                   ; restore
  466.                      POP    SI                   ;  register variables
  467.                      POP    DS                   ; restore caller's data segm.
  468.                      POP    BP                   ; save caller's base pointer
  469.                      RET    4                    ; pop parameters and return
  470.          APPLE87     ENDP
  471.  
  472.          CODE        ENDS
  473.  
  474.                      END
  475.  
  476.          ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  477.  
  478.          UNIT Time;
  479.  
  480.          INTERFACE
  481.  
  482.          FUNCTION Clock: LONGINT;          { same as VMS; time in milliseconds }
  483.  
  484.  
  485.          IMPLEMENTATION
  486.  
  487.          FUNCTION Clock: LONGINT; ASSEMBLER;
  488.          ASM
  489.                     PUSH    DS            { save caller's data segment }
  490.                     XOR     DX, DX        { initialize data segment to }
  491.                     MOV     DS, DX        {  access ticker counter }
  492.                     MOV     BX, 46Ch      { offset of ticker counter in segm.}
  493.                     MOV     DX, 43h       { timer chip control port }
  494.                     MOV     AL, 4         { freeze timer 0 }
  495.                     PUSHF                 { save caller's int flag setting }
  496.                     STI                   { allow update of ticker counter }
  497.                     LES     DI, DS:[BX]   { read BIOS ticker counter }
  498.                     OUT     DX, AL        { latch timer 0 }
  499.                     LDS     SI, DS:[BX]   { read BIOS ticker counter }
  500.                     IN      AL, 40h       { read latched timer 0 lo-byte }
  501.                     MOV     AH, AL        { save lo-byte }
  502.                     IN      AL, 40h       { read latched timer 0 hi-byte }
  503.                     POPF                  { restore caller's int flag }
  504.                     XCHG    AL, AH        { correct order of hi and lo }
  505.                     MOV     CX, ES        { ticker counter 1 in CX:DI:AX }
  506.                     CMP     DI, SI        { ticker counter updated ? }
  507.                     JE      @no_update    { no }
  508.                     OR      AX, AX        { update before timer freeze ? }
  509.                     JNS     @no_update    { no }
  510.                     MOV     DI, SI        { use second }
  511.                     MOV     CX, DS        {  ticker counter }
  512.          @no_update:NOT     AX            { counter counts down }
  513.                     MOV     BX, 36EDh     { load multiplier }
  514.                     MUL     BX            { W1 * M }
  515.                     MOV     SI, DX        { save W1 * M (hi) }
  516.                     MOV     AX, BX        { get M }
  517.                     MUL     DI            { W2 * M }
  518.                     XCHG    BX, AX        { AX = M, BX = W2 * M (lo) }
  519.                     MOV     DI, DX        { DI = W2 * M (hi) }
  520.                     ADD     BX, SI        { accumulate }
  521.                     ADC     DI, 0         {  result }
  522.                     XOR     SI, SI        { load zero }
  523.                     MUL     CX            { W3 * M }
  524.                     ADD     AX, DI        { accumulate }
  525.                     ADC     DX, SI        {  result in DX:AX:BX }
  526.                     MOV     DH, DL        { move result }
  527.                     MOV     DL, AH        {  from DL:AX:BX }
  528.                     MOV     AH, AL        {   to }
  529.                     MOV     AL, BH        {    DX:AX:BH }
  530.                     MOV     DI, DX        { save result }
  531.                     MOV     CX, AX        {  in DI:CX }
  532.                     MOV     AX, 25110     { calculate correction }
  533.                     MUL     DX            {  factor }
  534.                     SUB     CX, DX        { subtract correction }
  535.                     SBB     DI, SI        {  factor }
  536.                     XCHG    AX, CX        { result back }
  537.                     MOV     DX, DI        {  to DX:AX }
  538.                     POP     DS            { restore caller's data segment }
  539.          END;
  540.  
  541.  
  542.          BEGIN
  543.             Port [$43] := $34;           { need rate generator, not square wave}
  544.             Port [$40] := 0;             { generator as prog. by some BIOSes }
  545.             Port [$40] := 0;             { for timer 0 }
  546.          END. { Time }
  547.  
  548.  
  549.          ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  550.  
  551.          {$A+,B-,R-,I-,V-,N+,E+}
  552.          PROGRAM PeakFlop;
  553.  
  554.          USES Time;
  555.  
  556.          TYPE ParamRec = RECORD
  557.                             MaxIter, MaxRad, YMax, XMax: WORD;
  558.                             Qmax, Qmin, Pmax, Pmin: DOUBLE;
  559.                             FastMod: WORD;
  560.                             PlotFkt: POINTER;
  561.                             FLOPS:LONGINT;
  562.                          END;
  563.  
  564.          VAR Param: ParamRec;
  565.              Start: LONGINT;
  566.  
  567.  
  568.          {$L APFELM4.OBJ}
  569.  
  570.          PROCEDURE Apple87 (VAR Param: ParamRec);     EXTERNAL;
  571.  
  572.  
  573.          BEGIN
  574.             WITH Param DO BEGIN
  575.                MaxIter:= 50;
  576.                MaxRad := 30;
  577.                YMax   := 30;
  578.                XMax   := 30;
  579.                Pmin   :=-2.1;
  580.                Pmax   := 1.1;
  581.                Qmin   :=-1.2;
  582.                Qmax   := 1.2;
  583.                FastMod:= Word (FALSE);
  584.                PlotFkt:= NIL;
  585.                Flops  := 0;
  586.             END;
  587.             Start := Clock;
  588.             Apple87 (Param);         { executes 104002 FLOP }
  589.             Start := Clock - Start;  { elapsed time in milliseconds }
  590.             WriteLn ('Peak-MFLOPS: ', 104.002 / Start);
  591.          END.
  592.  
  593.          ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  594.  
  595.          ; FILE: M4X4.ASM
  596.          ;
  597.          ; assemble with TASM /e M4X4 or MASM /e M4X4
  598.  
  599.          CODE      SEGMENT BYTE PUBLIC 'CODE'
  600.  
  601.                    ASSUME  CS:CODE
  602.  
  603.                    PUBLIC  MUL_4x4
  604.                    PUBLIC  IIT_MUL_4x4
  605.  
  606.  
  607.          FSBP0     EQU     DB  0DBh, 0E8h        ; declare special IIT
  608.          FSBP1     EQU     DB  0DBh, 0EBh        ;  instructions
  609.          FSBP2     EQU     DB  0DBh, 0EAh
  610.          F4X4      EQU     DB  0DBh, 0F1h
  611.  
  612.  
  613.          ;---------------------------------------------------------------------
  614.          ;
  615.          ; MUL_4x4 multiplicates a four-by-four matrix by an array of four
  616.          ; dimensional vectors. This operation is needed for 3D transformations
  617.          ; in graphics data processing. There are arrays for each component of
  618.          ; a vector. Thus there is an ; array containing all the x components,
  619.          ; another containing all the y components and so on. Each component is
  620.          ; an 8 byte IEEE floating point number. Two indices into the array of
  621.          ; vectors are given. The first is the index of the vector that will be
  622.          ; processed first, the second is the index of the vector processed
  623.          ; last.
  624.          ;
  625.          ;---------------------------------------------------------------------
  626.  
  627.          MUL_4x4   PROC    NEAR
  628.  
  629.                    AddrX   EQU DWORD PTR [BP+24] ; address of X component array
  630.                    AddrY   EQU DWORD PTR [BP+20] ; address of Y component array
  631.                    AddrZ   EQU DWORD PTR [BP+16] ; address of Z component array
  632.                    AddrW   EQU DWORD PTR [BP+12] ; address of W component array
  633.                    AddrT   EQU DWORD PTR [BP+8]  ; addr. of 4x4 transform. mat.
  634.                    F       EQU WORD  PTR [BP+6]  ; first vector to process
  635.                    K       EQU WORD  PTR [BP+4]  ; last vector to process
  636.                    RetAddr EQU WORD  PTR [BP+2]  ; return address saved by call
  637.                    SavdBP  EQU WORD  PTR [BP+0]  ; saved frame pointer
  638.                    SavdDS  EQU WORD  PTR [BP-2]  ; caller's data segment
  639.  
  640.                    PUSH    BP                    ; save TURBO-Pascal frame pointer
  641.                    MOV     BP, SP                ; new frame pointer
  642.                    PUSH    DS                    ; save TURBO-Pascal data segment
  643.  
  644.                    MOV     CX, K                 ; final index
  645.                    SUB     CX, F                 ; final index - start index
  646.                    JNC     $ok                   ; must not
  647.                    JMP     $nothing              ;  be negative
  648.          $ok:      INC     CX                    ; number of elements
  649.  
  650.                    MOV     SI, F                 ; init offset into arrays
  651.                    SHL     SI, 1                 ; each
  652.                    SHL     SI, 1                 ;  element
  653.                    SHL     SI, 1                 ;   has 8 bytes
  654.  
  655.                    LDS     DI, AddrT             ; addr. of transformation mat.
  656.                    FLD     QWORD PTR [DI]        ; load a[0,0]   = R7
  657.                    FLD     QWORD PTR [DI+8]      ; load a[0,1]   = R6
  658.  
  659.          $mat_mul: LES     BX, AddrX             ; addr. of x component array
  660.                    FLD     QWORD PTR ES:[BX+SI]  ; load x[a]     = R5
  661.                    LES     BX, AddrY             ; addr. of y component array
  662.                    FLD     QWORD PTR ES:[BX+SI]  ; load y[a]     = R4
  663.                    LES     BX, AddrZ             ; addr. of z component array
  664.                    FLD     QWORD PTR ES:[BX+SI]  ; load z[a]     = R3
  665.                    LES     BX, AddrW             ; addr. of w component array
  666.                    FLD     QWORD PTR ES:[BX+SI]  ; load w[a]     = R2
  667.  
  668.                    FLD     ST(5)                 ; load a[0,0]   = R1
  669.                    FMUL    ST, ST(4)             ; a[0,0] * x[a] = R1
  670.                    FLD     ST(5)                 ; load a[0,1]   = R0
  671.                    FMUL    ST, ST(4)             ; a[0,1] * y[a] = R0
  672.                    FADDP   ST(1), ST             ; a[0,0]*x[a]+a[0,1]*y[a]=R1
  673.                    FLD     QWORD PTR [DI+16]     ; load a[0,2]   = R0
  674.                    FMUL    ST, ST(3)             ; a[0,2] * z[a] = R0
  675.                    FADDP   ST(1), ST             ; a[0,0]*x[a]...a[0,2]*z[a]=R1
  676.                    FLD     QWORD PTR [DI+24]     ; load a[0,3]   = R0
  677.                    FMUL    ST, ST(2)             ; a[0,3] * w[a] = R0
  678.                    FADDP   ST(1), ST             ; a[0,0]*x[a]...a[0,3]*w[a]=R1
  679.                    LES     BX, AddrX             ; get address of x vector
  680.                    FSTP    QWORD PTR ES:[BX+SI]  ; write new x[a]
  681.  
  682.                    FLD     QWORD PTR [DI+32]     ; load a[1,0]   = R1
  683.                    FMUL    ST, ST(4)             ; a[1,0] * x[a] = R1
  684.                    FLD     QWORD PTR [DI+40]     ; load a[1,1]   = R0
  685.                    FMUL    ST, ST(4)             ; a[1,1] * y[a] = R0
  686.                    FADDP   ST(1), ST             ; a[1,0]*x[a]+a[1,1]*y[a]=R1
  687.                    FLD     QWORD PTR [DI+48]     ; load a[1,2]   = R0
  688.                    FMUL    ST, ST(3)             ; a[1,2] * z[a] = R0
  689.                    FADDP   ST(1), ST             ; a[1,0]*x[a]...a[1,2]*z[a]=R1
  690.                    FLD     QWORD PTR [DI+56]     ; load a[1,3]   = R0
  691.                    FMUL    ST, ST(2)             ; a[1,3] * w[a] = R0
  692.                    FADDP   ST(1), ST             ; a[1,0]*x[a]...a[1,3]*w[a]=R1
  693.                    LES     BX, AddrY             ; get address of y vector
  694.                    FSTP    QWORD PTR ES:[BX+SI]  ; write new y[a]
  695.  
  696.                    FLD     QWORD PTR [DI+64]     ; load a[2,0]   = R1
  697.                    FMUL    ST, ST(4)             ; a[2,0] * x[a] = R1
  698.                    FLD     QWORD PTR [DI+72]     ; load a[2,1]   = R0
  699.                    FMUL    ST, ST(4)             ; a[2,1] * y[a] = R0
  700.                    FADDP   ST(1), ST             ; a[2,0]*x[a]+a[2,1]*y[a]=R1
  701.                    FLD     QWORD PTR [DI+80]     ; load a[2,2]   = R0
  702.                    FMUL    ST, ST(3)             ; a[2,2] * z[a] = R0
  703.                    FADDP   ST(1), ST             ; a[2,0]*x[a]...a[2,2]*z[a]=R1
  704.                    FLD     QWORD PTR [DI+88]     ; load a[2,3]   = R0
  705.                    FMUL    ST, ST(2)             ; a[2,3] * w[a] = R0
  706.                    FADDP   ST(1), ST             ; a[2,0]*x[a]...a[2,3]*w[a]=R1
  707.                    LES     BX, AddrZ             ; get address of z vector
  708.                    FSTP    QWORD PTR ES:[BX+SI]  ; write new z[a]
  709.  
  710.                    FLD     QWORD PTR [DI+96]     ; load a[3,0]   = R1
  711.                    FMULP   ST(4), ST             ; a[3,0] * x[a] = R5
  712.                    FLD     QWORD PTR [DI+104]    ; load a[3,1]   = R1
  713.                    FMULP   ST(3), ST             ; a[3,1] * y[a] = R4
  714.                    FLD     QWORD PTR [DI+112]    ; load a[3,2]   = R1
  715.                    FMULP   ST(2), ST             ; a[3,2] * z[a] = R3
  716.                    FLD     QWORD PTR [DI+120]    ; load a[3,3]   = R1
  717.                    FMULP   ST(1), ST             ; a[3,3] * w[a] = R2
  718.                    FADDP   ST(1), ST             ; a[3,3]*w[a]+a[3,2]*z[a]=R3
  719.                    FADDP   ST(1), ST             ; a[3,3]*w[a]...a[3,1]*y[a]=R4
  720.                    FADDP   ST(1), ST             ; a[3,3]*w[a]...a[3,0]*x[a]=R5
  721.                    LES     BX, AddrW             ; get address of w vector
  722.                    FSTP    QWORD PTR ES:[BX+SI]  ; write new w[a]
  723.  
  724.                    ADD     SI, 8                 ; new offset into arrays
  725.                    DEC     CX                    ; decrement element counter
  726.                    JZ      $done                 ; no elements left, done
  727.                    JMP     $mat_mul              ; transform next vector
  728.  
  729.          $done:    FSTP     ST(0)                ; clear
  730.                    FSTP     ST(0)                ;  FPU stack
  731.          $nothing: POP      DS                   ; restore TP data segment
  732.                    POP      BP                   ; restore TP frame pointer
  733.                    RET      24                   ; pop parameters and return
  734.  
  735.          MUL_4X4   ENDP
  736.  
  737.  
  738.          ;---------------------------------------------------------------------
  739.          ;
  740.          ; IIT_MUL_4x4 multiplicates a four-by-four matrix by an array of four
  741.          ; dimensional vectors. This operation is needed for 3D transformations
  742.          ; in graphics data processing. There are arrays for each component of
  743.          ; a vector.  Thus there is an array containing all the x components,
  744.          ; another containing all the y components and so on. Each component is
  745.          ; an 8 byte IEEE floating point number. Two indices into the array of
  746.          ; vectors are given. The first is the index of the vector that will be
  747.          ; processed first, the second is the index of the vector processed
  748.          ; last. This subroutine uses the special instructions only available
  749.          ; on IIT coprocessors to provide fast matrix multiply capabilities.
  750.          ; So make sure to use it only on IIT coprocessors.
  751.          ;
  752.          ;---------------------------------------------------------------------
  753.  
  754.          IIT_MUL_4x4   PROC    NEAR
  755.  
  756.                    AddrX   EQU DWORD PTR [BP+24] ; address of X component array
  757.                    AddrY   EQU DWORD PTR [BP+20] ; address of Y component array
  758.                    AddrZ   EQU DWORD PTR [BP+16] ; address of Z component array
  759.                    AddrW   EQU DWORD PTR [BP+12] ; address of W component array
  760.                    AddrT   EQU DWORD PTR [BP+8]  ; addr. of 4x4 transf. matrix
  761.                    F       EQU WORD  PTR [BP+6]  ; first vector to process
  762.                    K       EQU WORD  PTR [BP+4]  ; last vector to process
  763.                    RetAddr EQU WORD  PTR [BP+2]  ; return address saved by call
  764.                    SavdBP  EQU WORD  PTR [BP+0]  ; saved frame pointer
  765.                    SavdDS  EQU WORD  PTR [BP-2]  ; caller's data segment
  766.                    Ctrl87  EQU WORD  PTR [BP-4]  ; caller's 80x87 control word
  767.  
  768.                    PUSH    BP                    ; save TURBO-Pascal frame ptr
  769.                    MOV     BP, SP                ; new frame pointer
  770.                    PUSH    DS                    ; save TURBO-Pascal data seg.
  771.                    SUB     SP, 2                 ; make local variabe
  772.                    FSTCW   [Ctrl87]              ; save 80x87 ctrl word
  773.                    LES     SI, AddrT             ; ptr to transformation matrix
  774.                    FINIT                         ; initialize coprocessor
  775.                    FSBP2                         ; set register bank 2
  776.                    FLD     QWORD PTR ES:[SI]     ; load a[0,0]
  777.                    FLD     QWORD PTR ES:[SI+32]  ; load a[1,0]
  778.                    FLD     QWORD PTR ES:[SI+64]  ; load a[2,0]
  779.                    FLD     QWORD PTR ES:[SI+96]  ; load a[3,0]
  780.                    FLD     QWORD PTR ES:[SI+8]   ; load a[0,1]
  781.                    FLD     QWORD PTR ES:[SI+40]  ; load a[1,1]
  782.                    FLD     QWORD PTR ES:[SI+72]  ; load a[2,1]
  783.                    FLD     QWORD PTR ES:[SI+104] ; load a[3,1]
  784.                    FINIT                         ; initialize coprocessor
  785.                    FSBP1                         ; set register bank 1
  786.                    FLD     QWORD PTR ES:[SI+16]  ; load a[0,2]
  787.                    FLD     QWORD PTR ES:[SI+48]  ; load a[1,2]
  788.                    FLD     QWORD PTR ES:[SI+80]  ; load a[2,2]
  789.                    FLD     QWORD PTR ES:[SI+112] ; load a[3,2]
  790.                    FLD     QWORD PTR ES:[SI+24]  ; load a[0,3]
  791.                    FLD     QWORD PTR ES:[SI+56]  ; load a[1,3]
  792.                    FLD     QWORD PTR ES:[SI+88]  ; load a[2,3]
  793.                    FLD     QWORD PTR ES:[SI+120] ; load a[3,3]
  794.  
  795.                                         ; transformation matrix loaded
  796.  
  797.                    MOV     AX, F                 ; index of first vector
  798.                    MOV     DX, K                 ; index of last vector
  799.  
  800.                    MOV     BX, AX                ; index 1st vector to process
  801.                    MOV     CL, 3                 ; component has 8 (2**3) bytes
  802.                    SHL     BX, CL                ; compute offset into arrays
  803.  
  804.                    FINIT                         ; initialize coprocessor
  805.                    FSBP0                         ; set register bank 0
  806.  
  807.          $mat_loop:LES     SI, AddrW             ; addr. of W component array
  808.                    FLD     QWORD PTR ES:[SI+BX]  ; W component current vector
  809.                    LES     SI, AddrZ             ; addr. of Z component array
  810.                    FLD     QWORD PTR ES:[SI+BX]  ; Z component current vector
  811.                    LES     SI, AddrY             ; addr. of Y component array
  812.                    FLD     QWORD PTR ES:[SI+BX]  ; Y component current vector
  813.                    LES     SI, AddrX             ; addr. of X component array
  814.                    FLD     QWORD PTR ES:[SI+BX]  ; X component current vector
  815.                    F4X4                          ; mul 4x4 matrix by 4x1 vector
  816.                    INC     AX                    ; next vector
  817.                    MOV     DI, AX                ; next vector
  818.                    SHL     DI, CL                ; offset of vector into arrays
  819.  
  820.                    FSTP    QWORD PTR ES:[SI+BX]  ; store X comp. of curr. vect.
  821.                    LES     SI, AddrY             ; address of Y component array
  822.                    FSTP    QWORD PTR ES:[SI+BX]  ; store Y comp. of curr. vect.
  823.                    LES     SI, AddrZ             ; address of Z component array
  824.                    FSTP    QWORD PTR ES:[SI+BX]  ; store Z comp. of curr. vect.
  825.                    LES     SI, AddrW             ; address of W component array
  826.                    FSTP    QWORD PTR ES:[SI+BX]  ; store W comp. of curr. vect.
  827.  
  828.                    MOV     BX, DI                ; ofs nxt vect. in comp. arrays
  829.                    CMP     AX, DX                ; nxt vector past upper bound?
  830.                    JLE     $mat_loop             ; no, transform next vector
  831.                    FLDCW   [Ctrl87]              ; restore orig 80x87 ctrl word
  832.  
  833.                    ADD      SP, 2                ; get rid of local variable
  834.                    POP      DS                   ; restore TP data segment
  835.                    POP      BP                   ; restore TP frame pointer
  836.                    RET      24                   ; pop parameters and return
  837.          IIT_MUL_4x4   ENDP
  838.  
  839.          CODE      ENDS
  840.  
  841.                    END
  842.  
  843.          ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  844.  
  845.          {$N+,E+}
  846.  
  847.          PROGRAM Trnsform;
  848.  
  849.          USES Time;
  850.  
  851.          CONST VectorLen = 8190;
  852.  
  853.          TYPE  Vector    = ARRAY [0..VectorLen] OF DOUBLE;
  854.                VectorPtr = ^Vector;
  855.                Mat4      = ARRAY [1..4, 1..4] OF DOUBLE;
  856.  
  857.          VAR   X: VectorPtr;
  858.                Y: VectorPtr;
  859.                Z: VectorPtr;
  860.                W: VectorPtr;
  861.                T: Mat4;
  862.                K: INTEGER;
  863.                L: INTEGER;
  864.                First: INTEGER;
  865.                Last:  INTEGER;
  866.                Start: LONGINT;
  867.                Elapsed:LONGINT;
  868.  
  869.          PROCEDURE MUL_4X4     (X, Y, Z, W: VectorPtr;
  870.                                 VAR T: Mat4; First, Last: INTEGER); EXTERNAL;
  871.          PROCEDURE IIT_MUL_4X4 (X, Y, Z, W: VectorPtr;
  872.                                 VAR T: Mat4; First, Last: INTEGER); EXTERNAL;
  873.  
  874.          {$L M4X4.OBJ}
  875.  
  876.          BEGIN
  877.             WriteLn ('Test8087 = ', Test8087);
  878.             New (X);
  879.             New (Y);
  880.             New (Z);
  881.             New (W);
  882.             FOR L := 1 TO VectorLen DO BEGIN
  883.                X^ [L] := Random;
  884.                Y^ [L] := Random;
  885.                Z^ [L] := Random;
  886.                W^ [L] := Random;
  887.             END;
  888.             X^ [0] := 1;
  889.             Y^ [0] := 1;
  890.             Z^ [0] := 1;
  891.             W^ [0] := 1;
  892.             FOR K := 1 TO 4 DO BEGIN
  893.                FOR L := 1 TO 4 DO BEGIN
  894.                   T [K, L] := (K-1)*4 + L;
  895.                END;
  896.             END;
  897.             First := 0;
  898.             Last  := 8190;
  899.             Start := Clock;
  900.             MUL_4X4 (X, Y, Z, W, T, First, Last);
  901.             { IIT_MUL_4X4 (X, Y, Z, W, T, First, Last); }
  902.             Elapsed := Clock - Start;
  903.             WriteLn ('Number of vectors: ', Last-First+1);
  904.             WriteLn ('Time: ', Elapsed, ' ms');
  905.             WriteLn ('Equivalent to ', (28.0*(Last-First+1)/1e6)/
  906.                      (Elapsed*1e-3):0:4, ' MFLOPS');
  907.             WriteLn;
  908.             WriteLn ('Last vector:');
  909.             WriteLn;
  910.             WriteLn (X^[Last]);
  911.             WriteLn (Y^[Last]);
  912.             WriteLn (Z^[Last]);
  913.             WriteLn (W^[Last]);
  914.          END.
  915.