SEARCH FORPRM TV FORDBL DOUBLE PRECISION ROUTINES ,6(2031) ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SUBTTL REVISION HISTORY COMMENT \ ***** Begin Revision History ***** 1101 SWG 16-Aug-79 Cleanup FORDAR for V6 - take out F40 and KA conditionals remove following macros: FLADD, DFA, DAM, DFS, DSM FLMUL, DFM, DMM, FLDIV, DFD, DDM, all of which are KA only library routines. Replace DFN in IDF macro with DMOVN. remove KI conditionals because no longer necessary 1351 EDS 16-Mar-81 Q10-04786 Fix TWOSEG and RELOC problems. 1405 DAW 6-Apr-81 Extended addressing support. 1464 DAW 12-May-81 New error message format. 1673 CDM 9-Sept-81 Added new double precision routines. (DDIM, DINT, DNINT, DPROD, IDNINT) 1713 BL 15-Sep-81 Added 'DFL.3'. 1721 CDM 15-Sep-81 Edited DDIM to eliminate DDIM(-INF,+INF) overflow message. 1724 BL 17-Sep-81 Clean up error messages. 1746 CDM 25-Sep-81 Simple change to DDIM. ***** End Revision History ***** \ PRGEND TITLE DACOS ARC SINE AND ARC COSINE FUNCTIONS ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. JUNE 19, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DACOS EXTERN DACOS. DACOS=DACOS. PRGEND TITLE DASIN ARC SINE AND ARC COSINE FUNCTIONS ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. JUNE 19, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DASIN EXTERN DASIN. DASIN=DASIN. PRGEND TITLE DASIN. ARC SINE AND ARC COSINE FUNCTIONS ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. JUNE 19, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DASIN(X) AND DACOS(X) ARE CALCULATED AS FOLLOWS ; LET R(G) = (G*(RP1+G*(RP2+G*(RP3+G*(RP4+G*RP5)))))/(Q0+G*(Q1 ; +G*(Q2+G*(Q3+G*(Q4+G))))) ; (RP1,RP2,RP3,RP4,RP5,Q0,Q1,Q2,Q3,Q4, AND Q5 ARE GIVEN BELOW) ; LET S = Y + Y*R(G) ; FOR SITUATIONS 1. AND 3. BELOW, ; G = Y**2, WHERE Y = ABS(X) ; FOR SITUATIONS 2. AND 4., G = (1.0 - ABS(X))/2.0 ; AND Y = -2.0*SQRT(G) ; LET W = ABS(X) AND TERM THE RESULT T. ; THEN, CONSIDER THE FOUR SITUATIONS: ; 1. DASIN, FOR W <= 0.5. T = S, AND T IS NEGATED ; FOR NEGATIVE X ; 2. DASIN, FOR W > 0.5. T = PI/2 + S, AND T IS NEGATED ; FOR NEGATIVE X ; 3. DACOS, FOR W <= 0.5. T = PI/2 - S IF X < 0 ; = PI + S IF X > 0 ; 4. DACOS FOR W > 0.5. T = -S IF X < 0 ; = PI + S IF X > 0. ;THE RANGE OF DEFINITION FOR DASIN/DACOS IS ABS(X) <= 1.0. AND ERROR MESSAGES ; WILL RESULT FOR ARGUMENTS OUT OF THAT RANGE. DASIN/DACOS WILL BE SET ; TO + MACHINE INFINITY. ;REQUIRED (CALLED) ROUTINES: DSQRT ;REGISTERS T2, T3, T4, AND T5 WERE SAVED, USED, AND RESTORED ;THE ROUTINE HAS THE FOLLOWING CALLING SEQUENCE: ; XMOVEI L,ARG ; PUSHJ P,DASIN ; OR ; PUSHJ P,DACOS ;THE HIGH ORDER PART OF THE ANSWER IS RETURNED IN ACCUMULATOR T0 ;THE LOW ORDER PART OF THE ANSWER IS RETURNED IN ACCUMULATOR T1 SEARCH FORPRM TWOSEG 400000 SALL HELLO (DACOS,.) ;ENTRY TO DACOS ROUTINE HRRZI T0,2 ;SET FLAG TO TWO MOVEM T0,FLAG ;MOVE FLAG TO MEMORY JRST GETX ;GO TO GETX HELLO (DASIN,.) ;ENTRY TO DASIN ROUTINE SETZM FLAG ;SET FLAG TO ZEROES GETX: DMOVE T0,@(L) ;OBTAIN X JUMPGE T0,XPOS ;IF X IS NEGATIVE DMOVN T0,T0 ;Y = -X XPOS: CAMGE T0,HALF ;IF Y IS .LT. 1/2 JRST LE ;GO TO LE CAME T0,HALF ;IF Y IS .GT. 1/2 JRST GT ;GO TO GT JUMPE T1,LE ; GT: CAMGE T0,ONE ;IF Y IS .LT. ONE JRST ALG ;GO TO ALG CAMN T0,ONE ;IF Y IS .GT. ONE, GO TO ERR JUMPE T1,ALG ERR0: LERR (LIB,%,DASIN or DACOS; ABS(arg) > one; result = +infinity) HRLOI T0,377777 ;RESULT = SETO T1, ;+MACHINE INFINITY GOODBY (1) ;RETURN ALG: PUSH P,T2 ;SAVE ACCUMULATORS PUSH P,T3 PUSH P,T4 PUSH P,T5 MOVN T5,FLAG ;GET A COPY OF - FLAG HRRZI T4,2 ;SET T4 TO 2 ADD T5,T4 ;I = 2-FLAG HRLZI T2,201400 ;SET T2,T3 TO SETZ T3, ;ONE DFSB T2,T0 ;G = 1-Y FSC T2,-1 ;* 1/2 DMOVEM T2,TEMP ;SAVE A COPY OF G FUNCT DSQRT., ;Y = DSQRT(G) FSC T0,1 ;* 2 DMOVN T0,T0 ;Y = -Y MOVEM T5,I ;MOVE I TO MEMORY JRST GOTY ;GO TO GOTY LE: PUSH P,T2 ;SAVE ACCUMULATORS PUSH P,T3 PUSH P,T4 PUSH P,T5 MOVE T5,FLAG ;I = FLAG MOVEM T5,I ;MOVE I TO MEMORY CAMGE T0,TNM11H ;IF Y < TNM11H JRST GOTRES ;GO TO GOTRES DMOVE T2,T0 ;SAVE A COPY OF Y DFMP T2,T2 ;G = Y*Y GOTY: DMOVE T4,T2 ;SAVE A COPY OF G DFAD T4,Q4 ;XDEN = G + Q4 DFMP T4,T2 ;*G DFAD T4,Q3 ;+Q3 DFMP T4,T2 ;*G DFAD T4,Q2 ;+Q2 DFMP T4,T2 ;*G DFAD T4,Q1 ;+Q1 DFMP T4,T2 ;*G DFAD T4,Q0 ;+Q0 DMOVEM T2,TEMP ;SAVE A COPY OF G DFMP T2,RP5 ;XNUM = G*RP5 DFAD T2,RP4 ;+RP4 DFMP T2,TEMP ;*G DFAD T2,RP3 ;+RP3 DFMP T2,TEMP ;*G DFAD T2,RP2 ;+RP2 DFMP T2,TEMP ;*G DFAD T2,RP1 ;+RP1 DFMP T2,TEMP ;*G DFDV T2,T4 ;RESULT = XNUM/XDEN DFMP T2,T0 ;*Y DFAD T0,T2 ;+Y GOTRES: MOVE T5,I ;GET A COPY OF I SKIPE FLAG ;IF FLAG IS .NE. 0 JRST NE ;GO TO NE DFAD T0,A(T5) ;RESULT = RESULT + A(I) SKIPGE @(L) ;IF X IS NEGATIVE DMOVN T0,T0 ;RESULT = -RESULT POP P,T5 ;RESTORE ACCUMULATORS POP P,T4 POP P,T3 POP P,T2 GOODBY (1) ;RETURN NE: SKIPGE @(L) ;IF X IS NEGATIVE JRST NEGX ;GO TO NEGX DMOVN T0,T0 ;RESULT = -RESULT DFAD T0,A(T5) ;+A(I) POP P,T5 ;RESTORE ACCUMULATORS POP P,T4 POP P,T3 POP P,T2 GOODBY (1) ;RETURN NEGX: DFAD T0,B(T5) ;RESULT = RESULT + B(I) POP P,T5 ;RESTORE ACCUMULATORS POP P,T4 POP P,T3 POP P,T2 GOODBY (1) ;RETURN RP1: DOUBLE 572112065225,361372635617 ;-.27368494524164255994D+2 RP2: DOUBLE 206711524715,202552212153 ;.57208227877891731407D+2 RP3: DOUBLE 571302372325,226222540767 ;-.39688862997504877339D+2 RP4: DOUBLE 204504702731,073034447126 ;.10152522233806463645D+2 RP5: DOUBLE 577233210222,206360633365 ;-.69674573447350646411D+0 Q0: DOUBLE 567267447760,165074063632 ;-.16421096714498560795D+3 Q1: DOUBLE 211641111704,007534102703 ;.41714430248260412556D+3 Q2: DOUBLE 566202106100,352253670303 ;-.38186303361750149284D+3 Q3: DOUBLE 210455717445,226216157457 ;.15095270841030604719D+3 Q4: DOUBLE 572202642744,101474740606 ;-.23823859153670238830D+2 TNM11H: 134537657770 ;HIGH ORDER PART OF 1.D-11 A: DOUBLE 0,0 A1: DOUBLE 201622077325,021026430215 ;PI/2 B: DOUBLE 202622077325,021026430215 ;PI B1: DOUBLE 201622077325,021026430215 ;PI/2 ONE: 201400000000 ;1.0 HALF: 200400000000 ;1/2 RELOC ;DATA TEMP: DOUBLE 0,0 I: 0 FLAG: 0 RELOC PRGEND TITLE DATAN2 ARC TAN FUNCTION ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. MAY 9, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DATAN2 EXTERN DATN2. DATAN2=DATN2. PRGEND TITLE DATAN ARC TAN FUNCTION ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. MAY 9, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DATAN EXTERN DATAN. DATAN=DATAN. PRGEND TITLE DATAN. ARC TAN FUNCTION ; (DOUBLE PRECISION FLOATING POINT) ;COPYRIGHT (C) 1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DATAN(X) is computed as follows: ; ;If X < 0, compute DATAN(|X|) below, then DATAN(X) = -DATAN(-X). ; ;If X > 0, use the identity ; ; DATAN(X) = DATAN(XHI) + DATAN(Z) ; Z = (X - XHI) / (1 + X*XHI) ; ;where XHI is chosen so that |Z| <= tan(pi/32). ; ;XHI is chosen to be exactly representable as a double precision number, ;and so Z can be calculated without loss of significance. ; ;DATAN(XHI) is found by table lookup. It is stored as ATANHI + ATANLO to ;provide guard bits for the final addition to DATAN(Z). ; ;DATAN(Z) is evaluated by means of a polynomial approximation from Hart et al. ;(formula 4904). ; ;If X < tan(pi/32), DATAN(X) = DATAN(Z). ;If X > 1/tan(pi/32), DATAN(X) = pi/2 - DATAN(1/X). ; ;If tan(pi/32) < X < 1/tan(pi/32), an appropriate XHI is obtained by indexing ;into a table. The table tells which XHI to use for various ranges of X. ;The index into the table is formed from the low 3 exponent bits and the high ;3 fraction bits of X. SEARCH FORPRM TWOSEG 400000 SALL T6=6 ;MORE AC NAMES T7=7 T8=10 HELLO (DATAN,.) ;DATAN ENTRY PUSH P,T2 ;SAVE REGISTERS PUSH P,T3 PUSH P,T4 PUSH P,T5 PUSH P,T6 DMOVE T0,@(L) ;GET ARGUMENT X MOVEM T0,SGNFLG ;SAVE ARGUMENT SIGN FOR RESULT CAIGE T0,0 ;GET |X| DMOVN T0,T0 CAML T0,MAXX ;IS X LARGE? JRST LARGEX ;YES, GO USE ATAN(X) = PI/2 - ATAN(1/X) SETZ T2, ;T2 WILL GET OFFSET INTO XHI TABLES CAMG T0,MINX ;IS X SMALL ENOUGH THAT NO ARG REDUCTION IS REQUIRED? JRST CALC ;YES, JOIN CALCULATION BELOW MOVE T3,T0 ;GET HIGH WORD OF X LSHC T2,9 ;GET EXPONENT, SHIFT HIGH FRACTION BIT ;(ALWAYS 1) INTO SIGN BIT OF T3 ASHC T2,3 ;GET THREE FRACTION BITS, LEAVING ;THE 1 BEHIND HLRZ T2,OFFSET-2000+24(T2) ;GET OFFSET INTO XHI TABLES DMOVE T3,T0 ;GET A COPY OF X DFSB T0,XHI(T2) ;GET X-XHI DFMP T3,XHI(T2) ; X*XHI DFAD T3,ONE ; 1 + X*XHI DFDV T0,T3 ; (X-XHI) / (1 + X*XHI) ;Here T0 has the reduced argument Z with |Z| <= tan(pi/32). T2 has the ;index into the ATAN(XHI) tables ATANHI and ATANLO. SGNFLG is set negative ;if the result should be negated. CALC: DMOVE T3,T0 ;GET |Z| CAIGE T3,0 DMOVN T3,T3 CAMG T3,EPS ;IS IT SMALL ENOUGH THAT ATAN(Z) = Z? JRST SMALLX ;YES, BE FAST, AVOID UNDERFLOW DFMP T3,T3 ;GET Z**2 DMOVE T5,P06 ;GET P(Z**2) DFMP T5,T3 DFAD T5,P05 DFMP T5,T3 DFAD T5,P04 DFMP T5,T3 DFAD T5,P03 DFMP T5,T3 DFAD T5,P02 DFMP T5,T3 DFAD T5,P01 DFMP T3,T5 DFMP T3,T0 ; * Z DFAD T0,T3 ; + Z = ATAN(Z) SMALLX: DFAD T0,ATANLO(T2) ; + ATAN(XHI) LOW DFAD T0,ATANHI(T2) ; + ATAN(XHI) HI = DATAN(X) SKIPG SGNFLG ;ATTACH SIGN TO RESULT DMOVN T0,T0 RET: POP P,T6 ;RETURN POP P,T5 POP P,T4 POP P,T3 POP P,T2 POPJ P, LARGEX: DMOVE T3,T0 ;GET -1/X DMOVN T0,ONE DFDV T0,T3 MOVEI T2,PI2OFFS ;GET OFFSET OF PI/2 JRST CALC ;GO COMPUTE PI/2 + ATAN(-1/X) SUBTTL DATAN2 (Y,X) ;To compute DATAN2(Y,X), let U = |Y| and V = |X|, and compute DATAN(U/V). ;Then find DATAN2(Y,X) based on the signs of Y and X as follows: ; ; X Y DATAN2(Y/X) ; ; pos pos DATAN(U/V) ; pos neg -DATAN(U/V) ; neg pos -(DATAN(U/V) - pi) ; neg neg DATAN(U/V) - pi ; ;The add of -pi is combined with the add of DATAN(XHI) which is the last step ;of the DATAN algorithm. ; ;The reduced argument for DATAN is Z = (U/V - XHI) / (1 + U/V * XHI). ;This is rewritten as (U - V*XHI) / (V + U*XHI). To accurately calculate the ;numerator, find VHI and VLO with ; ; V = VHI + VLO ; VHI has at most 27 significant bits ; VLO has at most 35 significant bits ; ;and choose XHI with at most 13 significant bits. Then VHI*XHI and VLO*XHI can ;be exactly represented as double precision numbers, and the numerator is ; ; U - V*XHI = (U - VHI*XHI) - VLO*XHI HELLO (DATN2.,DATAN2) ;DATAN2 ENTRY PUSH P,T2 ;SAVE REGISTERS PUSH P,T3 PUSH P,T4 PUSH P,T5 PUSH P,T6 DMOVE T0,@0(L) ;GET Y DFDV T0,@1(L) ;GET Y/X JFCL EXCEP ;OVERFLOW AND UNDERFLOW CAN OCCUR DMOVEM T0,SGNFLG ;RESULT SHOULD BE MULTIPLIED BY SGN(Y/X) DMOVE T3,T0 ;GET |Y/X|, ATAN ARG CAIGE T3,0 DMOVN T3,T3 CAML T3,MAXX ;IS |Y/X| LARGE? JRST LARGE2 ;YES, GO USE ATAN(Y/X) = PI/2 - ATAN(X/Y) SETZ T2, ;GET OFFSET INTO ATAN TABLES CAMG T3,MINX ;IS |Y/X| SMALL ENOUGH TO USE POLYNOMIAL DIRECTLY? JRST [DMOVE T0,T3 ;YES, DO SO JRST SMALL2] LSHC T2,9 ;GET INDEX INTO OFFSET TABLES ASHC T2,3 HLRZ T2,OFFSET-2000+24(T2) ;GET INDEX INTO XHI TABLES PUSH P,T7 ;SAVE MORE REGISTERS PUSH P,T8 DMOVE T0,@0(L) ;GET |Y| = U CAIGE T0,0 DMOVN T0,0 DMOVEM T0,USAVE ;SAVE FOR LATER DMOVE T3,@1(L) ;GET |X| = V CAIGE T3,0 DMOVN T3,T3 MOVE T5,T3 ;GET A COPY OF V DMOVE T7,T3 ;GET ANOTHER SETZ T6, ;GET HIGH 27 BITS OF V = VHI DFSB T7,T5 ;GET LOW 35 BITS OF V = VLO DFMP T5,XHI(T2) ;GET V*XHI = VHI * XHI DFMP T7,XHI(T2) ; + VLO * XHI DFSB T0,T5 ;GET (U - VHI*XHI) DFSB T0,T7 ; - VLO*XHI DMOVE T5,USAVE ;GET U BACK DFMP T5,XHI(T2) ;GET U * XHI DFAD T3,T5 ;GET V + U*XHI DFDV T0,T3 ;GET (U - V*XHI) / (V + U*XHI) POP P,T8 ;RESTORE REGISTERS POP P,T7 SMALL2: SKIPGE @1(L) ;IF SECOND ARG (X) IS NEGATIVE ADDI T2,MPIOFFS ; ADD -PI TO RESULT JRST CALC ;GO GET ATAN AND RETURN LARGE2: DMOVN T0,@1(L) ;GET -X/Y DFDV T0,@0(L) MOVMM T0,SGNFLG ;SET SGNFLG POSITIVE MOVEI T2,PI2OFFS ;ADD PI/2 TO RESULT IF FIRST ARG (Y) IS POSITIVE SKIPGE @0(L) MOVEI T2,MPI2OFFS ;ADD -PI/2 TO RESULT IF FIRST ARG (Y) IS NEGATIVE JRST CALC ;GO COMPUTE +/- PI/2 + ATAN(-X/Y) EXCEP: SKIPN @1(L) ;CHECK FOR DIVIDE BY 0 JRST DIVCHK ;IF DIVIDE CHECK, GO CHECK FOR ATAN2(0,0) JUMPN T0,OVER ;IF OVERFLOW, RESULT IS +/- PI/2 SKIPL @1(L) ;IF UNDERFLOW, CHECK SECOND ARGUMENT LERR (LIB,%,,,RET) ;IF SECOND ARG (X) POSITIVE, RESULT UNDERFLOWS DMOVN T0,MPI ;ELSE RESULT IS PI WITH SIGN OF FIRST ARG JRST YSIGN ;GO ATTACH SIGN AND RETURN DIVCHK: SKIPE @0(L) ;CHECK FOR BOTH ARGS 0 JRST OVER ;ATAN2(NONZERO,0) IS SAME AS OVERFLOW LERR (LIB,%,) SETZB T0,T1 ;RETURN 0 JRST RET OVER: DMOVE T0,PI2 ;OVERFLOW, RESULT IS PI/2 WITH SIGN OF FIRST ARG YSIGN: SKIPGE @0(L) ;ATTACH SIGN OF FIRST ARGUMENT (Y) DMOVN T0,T0 JRST RET ;RETURN SUBTTL TABLES ;This table is indexed by the low 3 exponent bits and the high 3 fraction bits ;of X, where MINX < X < MAXX. It gives the offset into XHI, ATANHI, and ATANLO ;of a suitable XHI. DEFINE OFFS (X) < IRP X, > RADIX 10 OFFSET: OFFS <1,1,1,2,2,3,3,4,4,4,5,5,5,6,6,7,7,7,8,8,8,9,9,9> OFFS <10,10,10,10,11,11,11,12,12,12,12,12,13,13,13,13> OFFS <13,13,14,14,14,14,14,14,14,14,14,14,14,14,14> RADIX 8 PI2OFFS=2*^D15 ;OFFSET INTO ATANHI AND ATANLO OF PI/2 MPIOFFS=2*^D16 ;OFFSET OF -PI MPI2OFFS=2*^D31 ;OFFSET OF -PI/2 XHI: EXP 000000000000,0 ; .0000000000000000000 not used EXP 175657740000,0 ; .1054534912109375000 X1 EXP 176407740000,0 ; .1288757324218750000 X2 EXP 176477740000,0 ; .1562194824218750000 EXP 176617640000,0 ; .1952209472656250000 EXP 176777400000,0 ; .2497558593750000000 EXP 177477540000,0 ; .3121948242187500000 EXP 177617200000,0 ; .3898925781250000000 EXP 177776300000,0 ; .4984130859375000000 EXP 200515740000,0 ; .6522216796875000000 EXP 200674040000,0 ; .8673095703125000000 EXP 201453440000,0 ; 1.170166015625000000 EXP 201645100000,0 ; 1.645019531250000000 EXP 202511040000,0 ; 2.570800781250000000 EXP 203527000000,0 ; 5.359375000000000000 X14 ATANHI: EXP 000000000000,000000000000 ; .0000000000000000000 ATAN(0) EXP 175656261520,251717717115 ; .1050651824695016827 ATAN(X1) EXP 176406373154,305103547431 ; .1281692623534198424 EXP 176475276500,226075230141 ; .1549669515088296698 EXP 176612661305,353341503017 ; .1927961221331547523 EXP 176765175625,123431412160 ; .2447488705187413410 EXP 177465675077,207755453362 ; .3026068193134274483 EXP 177574536624,346575570716 ; .3717628303965550497 EXP 177731362665,337061624571 ; .4623772720674932799 EXP 200447716236,341667321523 ; .5779354489017924548 EXP 200555632634,072362412013 ; .7144577222199199908 EXP 200672140433,153152445005 ; .8636495725719767213 EXP 201406227204,252435004161 ; 1.024591515455200379 EXP 201463117110,112456461235 ; 1.199822549393342831 EXP 201542714675,142761043345 ; 1.386328658261967261 ATAN(X14) PI2: EXP 201622077325,021026430215 ; 1.570796326794896619 PI/2 MPI: EXP 575155700452,356751347563 ;-3.141592653589793238 -PI EXP 575173246105,164207746145 ;-3.036527471120291556 -PI + ATAN(X1) EXP 575176220221,273215536144 ;-3.013423391236373396 -PI + ATAN(X2) EXP 575201554376,370255221171 ;-2.986625702080963569 EXP 575206433527,115527433724 ;-2.948796531456638486 EXP 575215150344,104133030272 ;-2.896843783071051897 EXP 575224470162,337747115121 ;-2.838985834276365790 EXP 575235354335,213631126655 ;-2.769829823193238188 EXP 575251036741,252657532242 ;-2.679215381522299958 EXP 575267664122,247327234107 ;-2.563657204688000784 EXP 575311247221,375446052165 ;-2.427134931369873248 EXP 575334330561,311604060764 ;-2.277943081017816517 EXP 575361014155,104167751653 ;-2.117001138134592860 EXP 576016720236,050401400603 ;-1.941770104196450407 EXP 576076516023,100703762713 ;-1.755263995327825977 -PI + ATAN(X14) EXP 576155700452,356751347563 ;-1.570796326794896619 -PI/2 ATANLO: EXP 000000000000,000000000000 ; .0000000000000000000 EXP 701136265552,021570267746 ;-.1105497383317563340E-19 EXP 075632563733,177450755317 ; .5435918950106886920E-20 EXP 077603757522,264676361123 ; .2053885961335437547E-19 EXP 077723405617,212207745036 ; .2474984160195000159E-19 EXP 700044202740,142243000710 ;-.2518569148307390217E-19 EXP 100702445617,347031410060 ; .4770635578227287401E-19 EXP 100414524241,267545144713 ; .2844597940229199763E-19 EXP 100625025353,262726507370 ; .4288548085099372680E-19 EXP 676255576133,135473517361 ;-.7162797697818810113E-19 EXP 676323642342,201773750626 ;-.6356616556133364840E-19 EXP 077413353720,052424757433 ; .1415925447541254247E-19 EXP 100744743210,367241225634 ; .5134543068800780929E-19 EXP 101423305250,056513174673 ; .5831512827058072309E-19 EXP 102617464410,310763545272 ; .1692382723892392993E-18 EXP 101611431424,134015603464 ; .8333742918520878328E-19 EXP 675166346353,243762174314 ;-.1666748583704175666E-18 EXP 102634261642,105051607712 ; .1746358738541957411E-18 EXP 103601511025,077765605721 ; .3266520381981663157E-18 EXP 676115710653,365124065054 ;-.9192589013278796940E-19 EXP 675060707135,225203171017 ;-.1961351253927427867E-18 EXP 675072766647,260206474405 ;-.1918605498534914687E-18 EXP 101716137637,073361174657 ; .9787193190895619423E-19 EXP 674134635612,010745612637 ;-.3550693134652264558E-18 EXP 700335536464,205477162116 ;-.1536916027087339636E-19 EXP 103746522614,251310022042 ; .4122184681426969927E-18 EXP 103760133656,162370070313 ; .4202802795595514454E-18 EXP 101457607713,122451564536 ; .6432483060209586270E-19 EXP 103567657506,360715220731 ; .3183514413117920163E-18 EXP 676000222177,166457565523 ;-.1083597300998368435E-18 EXP 074603276433,074574160521 ; .2563414018821732668E-20 EXP 676166346353,243762174314 ;-.8333742918520878328E-19 ;COEFFICIENTS OF APPROXIMATION POLYNOMIAL ATAN(X) = X*P(X**2) P06: EXP 175461762413,330264551705 ; .0747006049800000000 P05: EXP 602213603465,225464265500 ;-.0908796288218500000 P04: EXP 175707070366,207135266005 ; .1111109168530032000 P03: EXP 601333333333,310141050127 ;-.1428571421988482590 P02: EXP 176631463146,146202001672 ; .1999999999989370798 P01: EXP 600252525252,252525266207 ;-.3333333333333326895 ONE: EXP 201400000000,000000000000 ; 1.000000000000000000 EPS: EXP 142561156421 ;LARGEST X WITH DATAN(X) = X (HIGH WORD) ;(IE, ALL DOUBLE PRECISION X WITH HIGH ; WORD OF X <= EPS HAVE DATAN(X)=X). MINX: EXP 175623327343 ;TAN(PI/32) (HIGH WORD, ROUNDED DOWN) MAXX: EXP 204504715427 ;1/TAN(PI/32) (HIGH WORD, ROUNDED UP) RELOC SGNFLG: BLOCK 1 ;SIGN TO BE ATTACHED TO RESULT USAVE: BLOCK 2 ;TEMP STORAGE FOR U RELOC PRGEND TITLE DCOSH HYPERBOLIC COSINE FUNCTION ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. JUNE 7, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DCOSH EXTERN DCOSH. DCOSH=DCOSH. PRGEND TITLE DCOSH. HYPERBOLIC COSINE FUNCTION ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. JUNE 7, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DCOSH(X) IS CALCULATED AS FOLLOWS ; LET V BE APPROXIMATELY 2 SO THAT LNV AND ABS(X)+LN V ; CAN BE EXACTLY REPRESENTED WHEN X IS EXACTLY REPRESENTABLE. ; THEN, LETTING Y = ABS(X), FOR ; Y <= 88.029678, DCOSH = (EXP(Y)+EXP(-Y))/2, ; 88.029678 <= Y < 128 * LN(2) ; DCOSH = (V/2)*EXP(Y - LN V), ; Y >= 128 * LN(2), DCOSH = +MACHINE INFINITY ;THE RANGE OF DEFINITION FOR DCOSH IS ABS(X) < 128 * LN(2) AND ERROR MESSAGES ; WILL RESULT FOR ARGUMENTS OUT OF THAT RANGE. DCOSH WILL BE SET ; TO + MACHINE INFINITY. ;REQUIRED (CALLED) ROUTINES: DEXP ;REGISTERS T2, T3, T4, AND T5 WERE SAVED, USED, AND RESTORED ;THE ROUTINE HAS THE FOLLOWING CALLING SEQUENCE: ; XMOVEI L,ARG ; PUSHJ P,DCOSH ;THE HIGH ORDER PART OF THE ANSWER IS RETURNED IN ACCUMULATOR T0 ;THE LOW ORDER PART OF THE ANSWER IS RETURNED IN ACCUMULATOR T1 SEARCH FORPRM TWOSEG 400000 SALL HELLO (DCOSH,.) ;ENTRY FOR DCOSH ROUTINE DMOVE T0,@(L) ;OBTAIN X JUMPGE T0,DCSH ;IF X IS NEGATIVE DMOVN T0,T0 ;Y = -X DCSH: CAMGE T0,TWENT2 ;IF ABS(X) < 22, GO TO JRST ALG2 ; FULL CALCULATION CAMG T0,BIGX ;IF Y IS .LE. BIGX JRST EASY ;THEN GO TO EASY CAMGE T0,YMAXHI ;IF HI OF Y < YMAXHI JRST GETW ; GO TO GETW CAME T0,YMAXHI ;IF HI OF Y > YMAXHI JRST ERRD ;GO TO ERRD CAMGE T1,YMAXLO ;HI OF Y = YMAXHI JRST GETW ; IF LO OF Y .LE. YMAXLO, GO TO GETW ERRD: HRLOI T0,377777 ;SET RESULT TO SETO T1, ;MACHINE INFINITY LERR (LIB,%,DCOSH: result overflow) GOODBY (1) ;RETURN GETW: DFAD T0,LNV ;W = Y-LNV DMOVEM T0,TEMP ;MOVE W TO A TEMP REGISTER FUNCT DEXP., ;Z = EXP(W) DMOVEM T0,TEMP ;SAVE A COPY OF Z FMP T0,CON1 ;RESULT = Z*CON1 SETZ T1, ;ZERO SECOND WORD DFAD T0,TEMP ;+Z JFCL ERRD GOODBY (1) ;RETURN EASY: DMOVEM T0,TEMP ;MOVE Y TO TEMP FUNCT DEXP., ;GET ITS EXPONENT FSC T0,-1 ;DIV BY 2 GOODBY (1) ;RETURN ALG2: CAMG T0,TM32 ;IF ABS(X) < 2**(-32) JRST TINY ; THE RESULT = 1. PUSH P,T2 ;SAVE MORE ACCUMULATORS PUSH P,T3 PUSH P,T4 PUSH P,T5 DMOVEM T0,TEMP ;MOVE Y TO A TEMPORARY REGISTER FUNCT DEXP., ;Z = DEXP(Y) DMOVE T2,T0 ;SAVE A COPY OF Z HRLZI T4,201400 ;SET T4 TO 1.0 SETZ T5, ;CLEAR SECOND WORD DFDV T4,T2 ;1/Z DFAD T0,T4 ;+Z FSC T0,-1 ;*1/2 POP P,T5 POP P,T4 ;RESTORE ACCUMULATORS POP P,T3 POP P,T2 GOODBY (1) ;RETURN TINY: HRLZI T0,201400 ;RESULT IS 1. SETZ T1, GOODBY (1) ;RETURN ONE: 201400000000 ;1.0 BIGX: 207540074635 ;88.029691 YMAXHI: 207542710277 ;88.722839111672999604 YMAXLO: 276434757153 ; LNV: DOUBLE 577235067500,101343000000 ;-.693147180559947174D0 CON1: 120414520000 ;.186417721E-14 TWENT2: 205540000000 ;22 TM32: 141400000000 ;2**(-32) RELOC ;DATA TEMP: DOUBLE 0,0 RELOC PRGEND TITLE DEXP EXPONENTIAL FUNCTION ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL.INC APRIL 3, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DEXP EXTERN DEXP. DEXP=DEXP. PRGEND TITLE DEXP. EXPONENTIAL FUNCTION ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. APRIL 3, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DEXP(X) IS CALCULATED AS FOLLOWS ; IF X <= -89.415986292232944914, EXP = 0 ; IF X > 88.029691931113054295, EXP = +MACHINE INFINITY ; IF X = 0.0, EXP = 1 ; OTHERWISE, ; THE ARGUMENT REDUCTION IS: ; LET X1 = [X], THE GREATEST INTEGER IN X ; X2 = X - X1 ; N = THE NEAREST INTEGER TO X/LN(2) ; THE REDUCED ARGUMENT IS G = ((X1 - N*C1)+X2)-N*C2 ; WHERE C1 = .543 (OCTAL), ; AND C2 IS GIVEN BELOW ; THE CALCULATION IS: ; EXP = R(G)*2**(N+1) ; WHERE R(G) = 0.5 + G*P/(Q - G*P) ; P = ((((P2*G**2)+P1)*G**2)+P0)*G**2 ; Q = (((((Q3*G**2)+Q2)*G**2)+Q1)*G**2)+Q0 ; P0, P1, P2, Q0, Q1, Q2, AND Q3 ARE GIVEN BELOW AS ; XP0, XP1, XP2, XQ0, XQ1, XQ2, AND XQ3 . ;THE RANGE OF DEFINITION FOR DEXP IS GIVEN ABOVE, AND ERROR MESSAGES ; WILL RESULT FOR ARGUMENTS OUT OF THAT RANGE ;REQUIRED (CALLED) ROUTINES: NONE ;REGISTERS T2, T3, T4, AND T5 WERE SAVED, USED, AND RESTORED ;THE ROUTINE HAS THE FOLLOWING CALLING SEQUENCE: ; XMOVEI L,ARG ; PUSHJ P,DEXP ;THE HIGH ORDER PART OF THE ANSWER IS RETURNED IN ACCUMULATOR T0, ;THE LOW ORDER PART OF THE ANSWER IS RETURNED IN ACCUMULATOR T1. SEARCH FORPRM TWOSEG 400000 SALL HELLO (DEXP,.) ;ENTRY TO DEXP ROUTINE DMOVE T0,@(L) ;GET DP ARGUMENT CAMLE T0,SMLXHI ;IF HI OF X > SMLXHI JRST NXTCHK ; GO TO NXTCHK CAME T0,SMLXHI ;IF HI OF X < SMLXHI JRST OUT2 ; GO TO OUT2 CAMGE T1,SMLXLO ;HI PART = SMLXHI, JRST OUT2 ; IF LO OF X < SMLXLO, OUT2 NXTCHK: CAMGE T0,BIGXHI ;IF HI OF X < BIGXHI JRST EXP1 ; GO TO EXP1 CAME T0,BIGXHI ;IF HI OF X > BIGXHI JRST ERRD ; GO TO ERR CAMG T1,BIGXLO ;IF LO OF X .LE. BIGXLO, JRST EXP1 ; GO TO EXP1 ERRD: LERR (LIB,%,DEXP: result overflow) HRLOI T0,377777 ;DEXP = +MACHINE INFINITY HRLOI T1,377777 GOODBY (1) ;RETURN OUT2: LERR (LIB,%,DEXP: result underflow) MOVEI T0,0 ;EXP = 0 MOVEI T1,0 GOODBY (1) ;RETURN EXP1: JUMPE T0,[MOVSI T0,(1.0) ;RETURN 1.0 FOR ARG OF ZERO GOODBY (1)] ;EXIT PUSH P,T2 ;SAVE ACCUMULATORS PUSH P,T3 PUSH P,T4 PUSH P,T5 MOVM T2,T0 ;GET |ARGUMENT| IN T2 CAML T2,LN2OV2 ;IS IT .LT. (1/2)*LN(2)? JRST REDUCE ;IF SO, IT MUST BE REDUCED. SETZ T4, ;MAKE SCALE INDEX 0. MOVEM T4,SAVEN ; AND SAVE IT. DMOVE T4,T0 ;GET COPY OF ARGUMENT. JRST MERGE ;MERGE WITH MAIN FLOW. REDUCE: DMOVE T2,T0 ;GET COPY OF ARG DFMP T2,RNDLN2 ;ARG/LN2 FIXR T2,T2 ;NEAREST INTEGER = N FLTR T2,T2 ;FLOAT N MOVEI T3,0 ; FIX T4,T0 ;[X] FLTR T4,T4 ;X1 = FLOAT[X] MOVEI T5,0 DFSB T0,T4 ;X2 = X - X1 MOVEM T2,SAVEN ;SAVEN DFMP T2,C1 ;N*C1 DFAD T4,T2 ;X1 + (N*C1) DFAD T4,T0 ;+X2 MOVE T2,SAVEN ;RETRIEVE N MOVEI T3,0 ;ZERO SECOND WORD DFMP T2,C2 ;FORM N*C2 DFAD T4,T2 ;(N*C2)+X1+(N*C1)+X2 DMOVE T0,T4 ;SAVE G MOVM T2,T4 ;GET ABS(G) MERGE: CAML T2,TWOM32 ;IF REDUCED ARG IS>= 2**-32 JRST APPRX ; GO TO APPRX DFAD T0,ONE ;R(G) = 1. + G FSC T0,-1 ;*.5 JRST BRNCH ;GO TO BRNCH APPRX: DFMP T4,T4 ;Z = G*G DMOVE T2,T4 ;SAVE Z DFMP T4,XP2 ;Z*XP2 DFAD T4,XP1 ;+XP1 DFMP T4,T2 ;* Z DFAD T4,XP0 ;+ XP0 DFMP T0,T4 ;* G DMOVE T4,T2 ;SAVE Z DFMP T4,XQ3 ;XQ3*Z DFAD T4,XQ2 ;+XQ2 DFMP T4,T2 ;*Z DFAD T4,XQ1 ;+ XQ1 DFMP T4,T2 ;* Z DFAD T4,XQ0 ; + XQ0 DFSB T4,T0 ; XQ - G*XP DFDV T0,T4 ;(G*XP)/(XQ-G*XP) DFAD T0,XQ0 ; + .5 BRNCH: MOVE T4,SAVEN ;RETRIEVE N FIX T4,T4 ADDI T4,1 ;N = N+1 FSC T0,0(T4) ;ADD N TO THE EXPONENT POP P,T5 ;RESTORE ACCUMULATORS POP P,T4 POP P,T3 POP P,T2 RET: GOODBY (1) ; RETURN SMLXHI: 570232254036 ;-89.415986292232944914 SMLXLO: 301750634730 ; BIGXHI: 207540074636 ;88.029691931113054295 BIGXLO: 077042573256 ; LN2OV2: 177542710300 ;(1/2) * LN(2) RNDLN2: DOUBLE 201561250731,112701376057 ;1.44269504088896341 = 1/LN2 ONE: DOUBLE 201400000000,000000000000 ;1.0D0 TWOM32: DOUBLE 141400000000,000000000000 ;2**-32 C1: DOUBLE 577235000000,000000000000 ;-0.693359375D0 C2: DOUBLE 164675002027,030206250331 ;2.12194440054690583D-4 XP0: DOUBLE 177400000000,000000000000 ;0.250D0 XP1: DOUBLE 171760351374,212411245446 ;0.757531801594227767D-2 XP2: DOUBLE 162410550412,271511036101 ;0.315551927656846464D-4 XQ0: DOUBLE 200400000000,000000000000 ;0.5D0 XQ1: DOUBLE 174721345024,167754776545 ;0.568173026985512218D-1 XQ2: DOUBLE 166512741427,012152337316 ;0.631218943743985036D-3 XQ3: DOUBLE 154623154303,071202125117 ;0.751040283998700461D-6 RELOC ;DATA SAVEN: 0 RELOC PRGEND TITLE DEXP2. DOUBLE ** INTEGER EXPONENENTIATION SUBTTL CHRIS SMITH/CKS 28-Jan-80 ;COPYRIGHT (C) 1980,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM TWOSEG 400000 SALL HELLO (DEXP2,.) PUSH P,T2 ;SAVE TEMP REGISTERS PUSH P,T3 PUSH P,T4 DMOVE T0,[EXP 201400000000,0] ;SET RESULT TO 1 MOVM T4,@1(L) ;GET ABS(EXPONENT) JUMPE T4,D2RET ;EXPONENT ZERO, GO RETURN 1 DMOVE T2,@0(L) ;GET BASE JUMPN T2,D2POS ;BASE NOT ZERO, GO GET BASE**ABS(EXP) JRST D2ZERO ;BASE ZERO, GO HANDLE D2LP: DFMP T2,T2 ;SQUARE BASE JFCL D2OVUN ;CATCH OVERFLOW D2POS: TRNE T4,1 ;CHECK LOW BIT OF EXPONENT DFMP T0,T2 ;MULTIPLY ANSWER BY BASE JFCL D2OVUN ;CATCH OVERFLOW LSH T4,-1 ;DIVIDE EXPONENT BY 2 JUMPN T4,D2LP ;HANDLE ALL BITS OF EXPONENT SKIPL @1(L) ;NEGATIVE EXPONENT? JRST D2RET ;NO, WE HAVE RESULT DMOVE T2,[EXP 201400000000,0] ;YES, GET RECIPROCAL DFDV T2,T0 DMOVE T0,T2 D2RET: POP P,T4 ;RESTORE TEMP REGISTERS POP P,T3 POP P,T2 POPJ P, ;DONE D2ZERO: SKIPL T4,@1(L) ;BASE 0, RESULT IS 0 IF EXPONENT POSITIVE JRST D2RTZ ;POSITIVE EXPONENT, GO RETURN 0 JRST D2OVFL ;ELSE RESULT OVERFLOWS D2OVUN: SKIPG T4,@1(L) ;NEGATIVE EXPONENT? JRST D2UNFL ;YES, RESULT UNDERFLOWS D2OVFL: LERR (LIB,%,DEXP2: result overflow) HRLOI T0,377777 ;OVERFLOW, GUESS POSITIVE RESULT HRLOI T1,377777 SKIPL @0(L) ;CHECK SIGN OF BASE JRST D2RET ;BASE NONNEGATIVE, GO RETURN +INFINITY TRNE T4,1 ;ODD EXPONENT? DMOVN T0,T0 ;NEGATIVE**ODD, RETURN -INFINITY JRST D2RET D2UNFL: LERR (LIB,%,DEXP2: result underflow) D2RTZ: SETZB T0,T1 ;RETURN 0 JRST D2RET PRGEND TITLE DEXP3. POWER FUNCTION ; (DOUBLE PRECISION) SUBTTL IMSL, INC. JUNE 4, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. ;DEXP3 CALCULATES X**Y WHERE X AND Y ARE FLOATING POINT VALUES IN THE ;FOLLOWING RANGES ; 0.0 < X < + MACHINE INFINITY (X MAY EQUAL 0.0 IF Y > 0.0 AND ; X MAY BE LESS THAN 0.0 IF Y IS AN INTEGER. IF X IS < 0.0 ; AND Y IS NOT AN INTEGER, A WARNING ERROR IS ISSUED AND ; DABS(X)**Y IS CALCULATED.) ; -128.9375 < FLOAT(INT((Y*LOG2(X))*16))/16 < 127.0 ; X**Y IS CALCULATED AS 2**W WHERE W = Y*LOG2(X). LOG2(X) IS ; CALCULATED AS FOLLOWS; ; X = F*(2**M), 1/2 <= F < 1. PICK P SUCH THAT P IS AN ODD ; INTEGER < 16 AND LET A = 2**(-P/16). NOW X = ((2**M)*A)*(F/A) ; LOG2(X) = M+LOG2(A) + LOG2(F/A) OR ; LOG2(X) = M-(P/16) + LOG2(F/A) . ; LET U1 = M-(P/16) AND ; U2 = LOG2(F/A) = LOG2((1+S)/(1-S)). ; THEN LOG2(X) = U1 + U2. ; AND S = (F-A)/(F+A). A RATIONAL ; APPROXIMATION IS USED TO EVALUATE U2. U1 AND U2 ARE THEN ; USED TO DETERMINE W1 AND W2 WHERE W=W1+W2 ; AND W1 = FLOAT(INT(W*16.0))/16.0. FINALLY Z=X**Y=2**W ; IS RECONSTRUCTED AS Z = (2**W1) * (2**W2) WHERE ; W1 = M1-P1/16 AND 2**W2 IS EVALUATED FROM ANOTHER RATIONAL ; FUNCTION. ;THE RANGE OF DEFINITION FOR DEXP3 IS GIVEN ABOVE, AND ERROR MESSAGES ; WILL RESULT FOR ARGUMENTS OUT OF THAT RANGE ;REQUIRED (CALLED) ROUTINES: NONE ;REGISTERS T2, T3, T4, T5, G1, G2, G3, AND G4 WERE SAVED, USED, AND RESTORED ;THE ROUTINE HAS THE FOLLOWING CALLING SEQUENCE: ; XMOVEI L,ARG ; PUSHJ P,EXP3 ;THE HIGH ORDER PART OF THE ANSWER IS RETURNED IN ACCUMULATOR T0 ;THE LOW ORDER PART OF THE ANSWER IS RETURNED IN T1 SEARCH FORPRM TWOSEG 400000 SALL HELLO (DEXP3,.) ;ENTRY TO DEXP3. ROUTINE PUSH P,T2 ;SAVE ACCUMULATORS PUSH P,T3 DMOVE T0,@(L) ;GET THE BASE DMOVE T2,@1(L) ;GET THE EXPONENT SETZM IFLAG ;SET INTEGER EXP FLAG TO 0 JUMPG T0,XOK ;IF X IS NOT GT 0 JUMPE T0,X0 ;AND IF X IS NOT = 0 PUSH P,T4 ;SAVE ACCUMULATORS PUSH P,T5 PUSH P,G1 PUSH P,G2 DMOVE T4,T2 ;GET A COPY OF Y JUMPGE T4,YP ;IF Y IS NEGATIVE DMOVN T4,T4 ;NEGATE IT YP: MOVE G2,T4 ;GET HIGH ORDER OF Y SETZ G1, ;SET G1 TO ZERO LSHC G1,11 ;GET EXPONENT OF Y SUBI G1,170 ;GET SHIFTING FACTOR LSHC T4,(G1) ;SHIFT OFF EXP AND PART OF INTEGER TLNE T4,400000 ;IF Y IS ODD SETOM IFLAG ;SET IFLAG TO ONES LSHC T4,1 ;SHIFT OFF LAST OF INTEGER PART DMOVN T0,T0 ;NEGATE X JUMPN T5,ERR0 ;IF LO OF Y .NE. 0, GO TO ERR JUMPE T4,YINT ;IF HI OF Y = 0, GO TO YINT ERR0: LERR (LIB,%,) JRST CONT ;GO TO CONT X0: JUMPG T2,RET1 ;IF Y IS NOT GT 0 LERR (LIB,%,) HRLOI T0,377777 ;RESULT = INFINITY HRLOI T1,377777 POP P,T3 ;RESTORE ACCUMULATORS POP P,T2 GOODBY (1) ;RETURN XOK: PUSH P,T4 ;SAVE ACCUMULATORS PUSH P,T5 PUSH P,G1 PUSH P,G2 YINT: JUMPN T3,CONT ;IF LO OF Y IS NOT 0, GO TO CONT JUMPN T2,YONE ;IF Y .NE. 0 GO TO YONE, OTHERWISE DMOVE T0,A1 ;SET RESULT TO 1.0 JRST RET2 ;GO TO RET2 YONE: CAMN T2,A1 ;IF Y = 1.0 JRST RET2 ;GO TO RET2 CONT: PUSH P,G3 PUSH P,G4 MOVE T4,T0 ;OBTAIN THE EXPONENT ASH T4,-33 ;SHIFT MANTISSA OFF SUBI T4,200 ;SUBTRACT 128 FROM EXPONENT AND T0,MASK1 ;EXTRACT MANTISSA IOR T0,MASK2 ;SET EXPONENT TO 0 ;THE FOLLOWING TESTS DETERMINE P MOVEI T5,2 ;NP = 1 CAMLE T0,A1+20 ;IF F GT A1(9) JRST NXT2 ;THEN GO TO NXT2 CAME T0,A1+20 ; JRST OK1 ; CAMLE T1,A1+21 ; JRST NXT2 OK1: ADDI T5,20 ; NXT2: CAMLE T0,A1+6(T5) JRST NXT4 ; CAME T0,A1+6(T5) ; JRST OK2 ; CAMLE T1,A1+7(T5) ; JRST NXT4 ; OK2: ADDI T5,10 ; NXT4: CAMLE T0,A1+2(T5) ; JRST NXT3 ; CAME T0,A1+2(T5) JRST OK3 ; CAMLE T1,A1+3(T5) ; JRST NXT3 ; OK3: ADDI T5,4 ; NXT3: MOVEM T5,PTEMP ;SAVE A COPY OF P DMOVE G1,T0 ;SAVE A COPY OF F DFAD G1,A1(T5) ;F+A1(P+1) DFSB T0,A1(T5) ;F-A1(P+1) SUBI T5,2 ;FORM P + 1 ASH T5,-1 ;(P+1)/2 DFSB T0,A2(T5) ;Z=(F-A1(P+1))-A2((P+1)/2) DFDV T0,G1 ;/(F+A1(P+1)) FSC T0,1 ; DMOVE G1,T0 ;SAVE A COPY OF Z DFMP G1,G1 ;FORM Z**2 DMOVE G3,G1 ;SAVE A COPY OF Z**2 DFMP G3,RP4 ;R(Z)=RP4*Z**2 DFAD G3,RP3 ;+RP3 DFMP G3,G1 ;*Z**2 DFAD G3,RP2 ;+RP2 DFMP G3,G1 ;*Z**2 DFAD G3,RP1 ;+RP1 DFMP G3,G1 ;*Z**2 DFMP G3,T0 ;*Z DFAD G3,T0 ;+Z DMOVE G1,G3 ;SAVE A COPY OF R(Z) DFMP G3,C ;U2 = R(Z) * C DFAD G3,G1 ;+ R(Z) ASH T4,4 ;U1 = M*16 MOVE T5,PTEMP ;GET A COPY OF P ASH T5,-1 SUB T4,T5 ;-P FLTR T4,T4 ;FLOAT IT SETZ T5, ; FSC T4,-4 ; DMOVE T0,T2 ;SAVE A COPY OF Y JUMPGE T2,YPOS ;IF Y IS NEGATIVE DMOVN T0,T0 ;NEGATE IT YPOS: MOVE G1,T0 ;GET COPY OF HIGH ORDER OF Y ASH G1,-33 ;SHIFT OFF MANTISSA CAIGE G1,227 ;IF THE EXPONENT IS LESS THAN 227 JRST SMALLY ;GO TO SMALLY CAIE G1,227 ;IF THE EXPONENT IS > 227 JRST LARGEY ;GO TO LARGEY SETZ T1, ;ZERO SECOND WORD JRST SGNCHK ;GO TO SGNCHK SMALLY: SUBI G1,175 ;GET SHIFTING FACTOR SETZ T1, ;SET Y1 TO 0 JUMPGE G1,GETY1 ;IF G1 IS .GE. 0 SETZ T0, ;THEN GO TO GETY1 JRST GETY2 ; GETY1: AND T0,MASK(G1) ;GET Y1 JRST SGNCHK ;GO TO CHECK SIGN LARGEY: CAIL G1,272 ;IF EXPONENT IS .GE. 272 JRST SGNCHK ;GO TO SGNCHK SUBI G1,230 ;GET SHIFTING FACTOR AND T1,REDMSK(G1) ;GET LOW ORDER PART OF Y1 SGNCHK: JUMPGE T2,GETY2 ;IF Y IS NEGATIVE DMOVN T0,T0 ;NEGATE Y1 GETY2: DMOVE G1,T2 ;SAVE A COPY OF Y DFSB T2,T0 ;Y2 = Y-Y1 DFMP T2,T4 ;FORM Y2*U1 JFCL DFMP G1,G3 ;FORM Y*U2 JFCL DFAD T2,G1 ;W = Y*U2+Y2*U1 DFMP T4,T0 ;FORM U1*Y1 JFCL DMOVE T0,T2 ;RECONSTRUCT W DFAD T0,T4 JFCL CAMG T0,BIGW ;IF W IS NOT TOO BIG JRST WOK ;GO TO WOK LERR (LIB,%,DEXP3: result overflow) HRLOI T0,377777 ; RESULT = INFINITY HRLOI T1,377777 JRST RET3 ;GO TO RET3 WOK: CAML T0,SMALLW ;IF W IS NOT TOO SMALL JRST WOK2 ;THEN PROCEED LERR (LIB,%,DEXP3: result underflow) SETZ T0, ; RESULT = 0 SETZ T1, POP P,G4 ;RESTORE ACCUMULATORS POP P,G3 POP P,G2 POP P,G1 POP P,T5 POP P,T4 POP P,T3 POP P,T2 GOODBY (1) ;RETURN WOK2: SETZ G2, ;ZERO G2 MOVM G1,T2 ;GET HIGH ORDER OF W MOVE G3,G1 ; ASH G3,-33 ;SHIFT OFF MANTISSA SUBI G3,175 ;GET SHIFTING FACTOR JUMPGE G3,GETW1 ;IF G3 .GE. 0 ;THEN GO TO GETW1 SETZ G1, ;OTHERWISE, SET W1 = 0 JRST GETW2 ;GO TO GETW2 GETW1: AND G1,MASK(G3) ;GETW1 JUMPG T2,GETW2 ;IF W IS NEGATIVE DMOVN G1,G1 ;NEGATE W1 GETW2: DFSB T2,G1 ;W2 = W-W1 DFAD T4,G1 ;W = W1+U1*Y1 MOVM T0,T4 ;SAVE A COPY OF ABS(W) MOVE G1,T0 ; SETZ T1, ;ZERO T1 ASH G1,-33 ;SHIFT OFF MANTISSA SUBI G1,175 ;GET SHIFTING FACTOR JUMPGE G1,GTW1 ;IF G1 IS .GE. 0 ;THEN GO TO GTW1 SETZ T0, ;OTHERWISE, SET W1 TO 0 JRST GTW2 ;GO TO GET W2 GTW1: AND T0,MASK(G1) ;GET W1 JUMPGE T4,GTW2 ;IF W IS NEGATIVE DMOVN T0,T0 ;NEGATE W1 GTW2: DFSB T4,T0 ;FORM W-W1 DFAD T2,T4 ;W2 = W2+(W-W1) MOVM T4,T2 ;SAVE A COPY OF ABS(W2) MOVE G1,T4 ; SETZ T5, ;ZERO T5 ASH G1,-33 ;SHIFT OFF MANTISSA SUBI G1,175 ;GET SHIFTING FACTOR JUMPGE G1,GW ;IF G1 IS .GE. 0 ;THEN GO TO GW SETZ T4, ;OTHERWISE, SET W=0 JRST GW2 ;GO TO GW2 GW: AND T4,MASK(G1) ;GET W JUMPGE T2,GW2 ;IF W2 IS NEGATIVE DMOVN T4,T4 ;NEGATE W GW2: DFAD T0,T4 ;FORM W1 + W FSC T0,4 ;*16 FIX T0,T0 ;IW1 DFSB T2,T4 ;W1 = W2 - W JUMPLE T2,W2POS ;IF W2 .GT. 0 DFSB T2,SXTNTH ;W2 = W2-.0625 ADDI T0,1 ;IW1 = IW1+1 W2POS: MOVE T5,T0 ;SAVE A COPY OF IW1 JUMPGE T5,NPOS ADDI T5,17 NPOS: ASH T5,-4 ;M1 = IW1/16 JUMPL T0,INEG ;IF IW1 .GE. 0 ADDI T5,1 ;M1 = M1+1 INEG: MOVE T4,T5 ;SAVE A COPY OF M1 ASH T4,4 ;P1 = 16*M1 SUB T4,T0 ;-IW1 ASH T4,1 ; DMOVE T0,T2 ;SAVE A COPY OF W2 DFMP T2,Q7 ;Z = Q7*W2 DFAD T2,Q6 ;+Q6 DFMP T2,T0 ;*W2 DFAD T2,Q5 ;+Q5 DFMP T2,T0 ;*W2 DFAD T2,Q4 ;+Q4 DFMP T2,T0 ;*W2 DFAD T2,Q3 ;+Q3 DFMP T2,T0 ;*W2 DFAD T2,Q2 ;+Q2 DFMP T2,T0 ;*W2 DFAD T2,Q1 ;+Q1 DFMP T0,T2 ;*W2 DFMP T0,A1(T4) ;*A1(P1+1) DFAD T0,A1(T4) ;+A1(P1+1) FSC T0,(T5) ;ADD M1 TO THE EXP OF Z RET3: POP P,G4 ;RESTORE ACCUMULATORS POP P,G3 RET2: SKIPE IFLAG ;IF Y IS ODD NEGATIVE INTEGER DMOVN T0,T0 ;NEGATE RESULT POP P,G2 POP P,G1 POP P,T5 POP P,T4 RET1: POP P,T3 POP P,T2 GOODBY (1) ;RETURN SXTNTH: DOUBLE 175400000000,000000000000 ;.0625D0 BIGW: 127.0E0 ;UPPER BOUND FOR WW SMALLW: -128.9375E0 ;LOWER BOUND FOR W RP1: DOUBLE 175525252525,125252514430 ;.833333333333332114D-1 RP2: DOUBLE 172631463146,147404016032 ;.125000000005037992D-01 RP3: DOUBLE 170444444413,211706653423 ;.223214212859242590D-2 RP4: DOUBLE 165707437566,322576050305 ;.434457756721631196D-3 C: DOUBLE 177705243545,053405770274 ;.442695040888963407D0 Q7: DOUBLE 160764733570,344546177513 ;.149288526805956082D-4 Q6: DOUBLE 164502757270,016046163745 ;.154002904409897646D-3 Q5: DOUBLE 167535417606,065746637065 ;.133335413135857847D-02 Q4: DOUBLE 172473125333,204616630130 ;.961812905951724170D-02 Q3: DOUBLE 174706541065,300601154766 ;.555041086640855953D-1 Q2: DOUBLE 176753767577,340542316147 ;.240226506959095371D0 Q1: DOUBLE 200542710277,276434757056 ;.693147180559945296D0 A1: DOUBLE 201400000000,000000000000 ;A1(I), I=1,17 = A12: DOUBLE 200752225750,251110331413 ;2**((1-I)/16). THIS A13: DOUBLE 200725403067,076722207114 ;TABLE IS SEARCHED A14: DOUBLE 200701463367,141251234104 ;TO DETERMINE P. A15: DOUBLE 200656423746,126551655275 ; A16: DOUBLE 200634222140,250770220071 ; A17: DOUBLE 200612634520,212520333267 ; A18: DOUBLE 200572042434,372600606660 ; A19: DOUBLE 200552023631,237635714441 ; A110: DOUBLE 200532540767,122052051262 ; A111: DOUBLE 200513773265,115425047073 ; A112: DOUBLE 200475724623,004432042153 ; A113: DOUBLE 200460337602,214333425134 ; A114: DOUBLE 200443417233,235261070316 ; A115: DOUBLE 200427127017,037250572672 ; A116: DOUBLE 200413253033,076304417305 ; A117: DOUBLE 200400000000,000000000000 ; A2: DOUBLE 077756350307,256663004531 A22: DOUBLE 101406261124,022156463124 A23: DOUBLE 702043260343,345371331424 A24: DOUBLE 676250402176,331550372175 A25: DOUBLE 676111401243,043654640102 A26: DOUBLE 101640442174,047535641316 A27: DOUBLE 676017655543,202562654434 A28: DOUBLE 100613445334,141025236276 MASK1: 000777777777 ;MASK FOR MANTISSA MASK2: 200000000000 ;MASK FOR EXPONENT REDMSK: 600000000000 RDMSK2: 700000000000 RDMSK3: 740000000000 RDMSK4: 760000000000 RDMSK5: 770000000000 RDMSK6: 774000000000 RDMSK7: 776000000000 RDMSK8: 777000000000 MASK: 777400000000 ;MASKS FOR FINDING W1 MSK1: 777600000000 ; MSK2: 777700000000 ; MSK3: 777740000000 ; MSK4: 777760000000 ; MSK5: 777770000000 ; MSK6: 777774000000 ; MSK7: 777776000000 ; MSK8: 777777000000 ; MSK9: 777777400000 ; MSK10: 777777600000 ; MSK11: 777777700000 MSK12: 777777740000 MSK13: 777777760000 MSK14: 777777770000 MSK15: 777777774000 MSK16: 777777776000 MSK17: 777777777000 MSK18: 777777777400 MSK19: 777777777600 MSK20: 777777777700 MSK21: 777777777740 MSK22: 777777777760 MSK23: 777777777770 MSK24: 777777777774 MSK25: 777777777776 MSK26: 777777777777 RELOC ;DATA PTEMP: 0 IFLAG: 0 RELOC PRGEND TITLE DLOG10 LOG BASE 10 FUNCTION ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. APRIL 8, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DLOG10 EXTERN DLG10. DLOG10=DLG10. PRGEND TITLE DLOG NATURAL LOG FUNCTION ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. APRIL 8, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DLOG EXTERN DLOG. DLOG=DLOG. PRGEND TITLE DLOG. NATURAL LOG FUNCTION ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. APRIL 8, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DLOG10(X) AND DLOG(X) ARE CALCULATED AS FOLLOWS ; FOR X = 0 AN ERROR MESSAGE IS ISSUED AND -MACHINE INFINITY ; IS RETURNED AS THE RESULT. ; FOR X < 0 AN ERROR MESSAGE IS ISSUED, X IS SET T0 -X AND ; CALCULATION CONTINUES. ; FOR X > 0, X = F*2**F(M), 1/2 < F < 1 ; DEFINE G AND N SO THAT F = G*2(-N), 1/SQRT(2) <= G < SQRT(2). ; NOW ; DLOG(X) = (K*M-N) * DLOG(2) + DLOG(G) ; AND ; DLOG10(X) = DLOG10(E) * DLOG(X) = DLOG(X)/DLOG(10) ; ; DLOG(G) IS EVALUATED BY DEFINING S = (G-1)/(G+1) AND Z = 2*S ; AND THEN CALCULATING DLOG(G) = DLOG((1+Z/2)/(1-Z/2)) USING ; A MINIMAX RATIONAL APPROXIMATION. ;THE RANGE OF DEFINITION FOR DLOG/DLOG10 IS GIVEN ABOVE, AND ERROR MESSAGES ; WILL RESULT FOR ARGUMENTS OUT OF THAT RANGE ;REQUIRED (CALLED) ROUTINES: NONE ;REGISTERS T2, T3, T4, AND T5 WERE SAVED, USED, AND RESTORED ;THE ROUTINE HAS THE FOLLOWING CALLING SEQUENCE: ; XMOVEI L,ARG ; PUSHJ P,DLOG ; OR ; PUSHJ P,DLOG10 ;THE HIGH ORDER PART OF THE ANSWER IS RETURNED IN ACCUMULATOR T0. ;THE LOW ORDER PART OF THE ANSWER IS RETURNED IN ACCUMULATOR T1. SEARCH FORPRM TWOSEG 400000 SALL HELLO (DLG10.,DLOG10) ;ENTRY TO LOG BASE 10 ROUTINE. SETZM FLAG ;CLEAR FLAG FOR DLOG10 ENTRY JRST ALG ;GO TO ALGORITHM HELLO (DLOG,.) ;ENTRY TO NATURAL LOG ROUTINE SETOM FLAG ;SET FLAG FOR DLOG ENTRY ALG: DMOVE T0,@(L) ;GET ARG JUMPG T0,STRT ;IF ARG > ZERO GO TO STRT JUMPN T0,ARGN ;OTHERWISE, IF ARG .NE. 0 GO TO ARGN LERR (LIB,%,DLOG or DLOG10: zero arg; result = -infinity) HRLZI T0,400000 ;SET RESULT TO HRRZI T1,000001 ;LARGE NEGATIVE NUMBER GOODBY (1) ARGN: LERR (LIB,%,DLOG or DLOG10: negative arg; result = LOG(ABS(arg))) DMOVN T0,T0 ;ARG = -ARG STRT: PUSH P,T2 ;SAVE ACCUMULATORS PUSH P,T3 PUSH P,T4 PUSH P,T5 MOVE T2,T0 ;GET COPY OF HIGH ORDER PART OF ARG ASH T2,-33 ;SHIFT OFF THE MANTISSA SUBI T2,200 ;SUBTRACT 128 FROM THE EXP DMOVE T4,T0 ;GET COPY OF ARG AND T4,MASK1 ;EXTRACT MANTISSA IOR T4,MASK2 ;SET EXP TO 0 CAMLE T4,HI ;IS HIGH PART OF F > SQRT(.5) JRST FGT ;YES, GO TO FGT CAME T4,HI ;NO, IS F = HIGH OF SQRT(.5) JRST FLT ;NO, GO TO FLT CAMLE T5,LOW ;YES, IS LOW PART > LOW PART OF SQRT(.5) JRST FGT ;YES, GO TO FGT FLT: SUBI T2,1 ;N = N-1 FLTR T2,T2 MOVEM T2,N ; SAVE N DFAD T4,MHALF ;ZNUM = F-.5 DMOVE T2,T4 ;GET COPY OF ZNUM FSC T4,-1 ;ZDEN = ZNUM * .5 DFAD T4,HALF ; + .5 JRST EVALRZ FGT: FLTR T2,T2 MOVEM T2,N ; SAVE N DMOVE T2,T4 DFAD T2,B3 ;ZNUM = F - 1.0 FSC T4,-1 ;ZDEN = F*.5 DFAD T4,HALF ; + .5 EVALRZ: DFDV T2,T4 ;Z = ZNUM/ZDEN DMOVE T4,T2 ; DFMP T4,T4 ;W = Z*Z DMOVE T0,T4 ;SAVE COPY OF W DFAD T4,B2 ; FORM B(W). B(W) = W + B2 DFMP T4,T0 ; * W DFAD T4,B1 ; + B1 DFMP T4,T0 ; * W DFAD T4,B0 ; + B0 DMOVEM T0,SAVEW ; SAVE A COPY OF W DFMP T0,A2 ;FORM A(W). A(W)= A2*W DFAD T0,A1 ; + A1 DFMP T0,SAVEW ; * W DFAD T0,A0 ; + A0 DFDV T0,T4 ;R(Z) = A(W)/B(W) DFMP T0,SAVEW ; * W DFMP T0,T2 ; *Z JFCL DFAD T0,T2 ; + Z MOVE T2,N ;RETRIEVE N MOVEI T3,0 ;ZERO OUT SECOND WORD DFMP T2,C1 ;FORM N*C1 MOVE T4,N ;RETRIEVE N MOVEI T5,0 ;ZERO OUT SECOND WORD DFMP T4,C2 ;FORM N*C2 DFAD T0,T4 ;RESULT = N*C2 + R(Z) DFAD T0,T2 ; + N*C1 SKIPN FLAG DFMP T0,C3 ;IF DLOG10 ROUTINE, RESULT = RESULT*C3 POP P,T5 ;RESTORE ACCUMULATORS POP P,T4 POP P,T3 POP P,T2 GOODBY (1) C1: DOUBLE 200543000000,000000000000 ;.693359375D0 C2: DOUBLE 613102775750,347571527443 ;-2.12194440054690583D-4 C3: DOUBLE 177674557305,111562416145 ;.434294481903251828D0 A0: DOUBLE 570377400073,123045145570 ;-.641249434237455811D2 A1: DOUBLE 205406111210,005527754477 ;.163839435630215342D2 A2: DOUBLE 577153575223,052241327104 ;-.789561128874912573D0 B0: DOUBLE 565177200130,374467610432 ;-.769499321084948798D3 B1: DOUBLE 211470020376,205223175724 ;.312032220919245328D3 B2: DOUBLE 571342517755,046030761117 ;-.356679777390346462D2 B3: DOUBLE 576400000000,000000000000 ;-.1D1 HALF: DOUBLE 200400000000,000000000000 ;.5D0 MHALF: DOUBLE 577400000000,000000000000 ;-.5D0 HI: 200552023631 LOW: 237635714441 MASK1: 000777777777 MASK2: 200000000000 RELOC ;DATA N: 0 FLAG: 0 SAVEW: DOUBLE 000000000000,000000000000 RELOC PRGEND TITLE DMOD DOUBLE PRECISION SUBTTL MARY PAYNE /MHP/CKS 25-Jan-80 ;COPYRIGHT (C) 1980,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. NOSYM ENTRY DMOD EXTERN DMOD. DMOD=DMOD. PRGEND TITLE DMOD. DOUBLE PRECISION SUBTTL MARY PAYNE /MHP/CKS 25-Jan-80 ;COPYRIGHT (C) 1980,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM TWOSEG 400000 SALL T6=6 ;EXTRA TEMP REGISTERS T7=7 T8=10 T9=11 HELLO (DMOD,.) DMOVEM T2,SAVE2 ;SAVE TEMP REGISTERS DMOVEM T4,SAVE4 DMOVEM T6,SAVE6 DMOVEM T8,SAVE8 DMOVE T0,@0(L) ;GET ABS(A) IN T0-T1 CAIGE T0,0 DMOVN T0,T0 REPEAT: DMOVE T2,@1(L) ;GET ABS(B) IN T2-T3 CAIGE T2,0 DMOVN T2,T2 CAMG T0,T2 ;CHECK IF A .LE. B JRST FAST ;YES, GO BE FAST DMOVEM T2,B ;SAVE ABS(B) DMOVE T8,T0 ;SAVE ABS(A) IN T8-T9 ;BREAK B INTO HIGH, MIDDLE, AND LOW PARTS ;B IN T2/ XXXHHHHHHHHH MMMMMMMMMLLL MOVE T6,T2 ;GET T6/ XXXHHHHHHHHH 000000000000 SETZ T7, ; DFSB T2,T6 ; T2/ XXXMMMMMMMMM LLL000000000 MOVE T4,T2 ; T4/ XXXMMMMMMMMM 000000000000 SETZ T5, ; DFSB T2,T4 ; T2/ XXXLLL000000 000000000000 DFDV T8,B ;GET A/B JFCL OVFL ;CATCH OVERFLOW CAML T8,[244400000000] ;IS QUOTIENT OVER 2**35? JRST [AND T9,[777000000000] ;YES, TRUNCATE IT TO 35 BITS JRST BIG] ; TO KEEP PRODUCTS EXACT FIX T9,T8 ;UNDER 2**35, TRUNCATE TO NEAREST INTEGER MOVSI T8,276000 ;DFLOAT THE INTEGER DFAD T8,[EXP 0,0] BIG: DFMP T6,T8 ;GET A-B*[A/B] IN T0-T1 DFSB T0,T6 DFMP T4,T8 DFSB T0,T4 DFMP T2,T8 DFSB T0,T2 OVCONT: CAIGE T0,0 ;IF RESULT NEGATIVE, [A/B] WAS 1 TOO HIGH DFAD T0,B ;CORRECT THE RESULT CAMGE T0,B ;IF RESULT .LT. B, OK JRST MRET CAMG T0,B ;ELSE [A/B] WAS TOO LOW CAML T1,B+1 ; AND MUST SUBTRACT SOME MORE MULTIPLES OF B JRST REPEAT ;GO GET T0 MOD B IN T0 MRET: SKIPGE @0(L) ;GIVE RESULT SIGN OF A DMOVN T0,T0 MRET1: DMOVE T8,SAVE8 ;RESTORE REGISTERS DMOVE T6,SAVE6 DMOVE T4,SAVE4 DMOVE T2,SAVE2 POPJ P, ;DONE FAST: CAMN T0,T2 ;HERE IF A HIGH .LE. B HIGH CAMGE T1,T3 ;HIGH WORDS EQUAL, COMPARE LOW WORDS JRST MRET ;A .LT. B, GO RETURN DMOD = A DFSB T0,T2 ;ELSE RETURN DMOD = A - B JRST MRET OVFL: JUMPE T6,RTZ ;IF B IS ZERO, DMOD IS 0 FSC T0,-177 ;SCALE A TO AVOID OVERFLOW IN DIVIDE DMOVE T8,T0 ;GET A DFDV T8,B ;GET A/B AND T9,[777000000000] ;TRUNCATE TO 35 BITS SO PRODUCTS ARE EXACT DFMP T6,T8 ;GET A-B*[A/B] DFSB T0,T6 DFMP T4,T8 DFSB T0,T4 DFMP T2,T8 DFSB T0,T2 FSC T0,177 ;UNDO SCALING JRST OVCONT ;CONTINUE RTZ: SETZB T0,T1 ;A .EQ. B, RETURN DMOD = 0 JRST MRET1 RELOC ;DATA B: BLOCK 2 ;ABS(B) SAVE2: BLOCK 2 ;TEMP REGISTERS SAVE4: BLOCK 2 SAVE6: BLOCK 2 SAVE8: BLOCK 2 RELOC PRGEND TITLE DCOS SINE AND COSINE FUNCTIONS ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. MAY 1, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DCOS EXTERN DCOS. DCOS=DCOS. PRGEND TITLE DSIN SINE AND COSINE FUNCTIONS ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. MAY 1, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DSIN EXTERN DSIN. DSIN=DSIN. PRGEND TITLE DSIN. SINE AND COSINE FUNCTIONS ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. MAY 1, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;COS(X) AND SIN(X) ARE CALCULATED AS FOLLOWS ; NOTE THAT COS(X) = COS(-X) AND SIN(X) = -SIN(-X) ; LET ABS(X) = N*PI + F WHERE ABS(F) <= PI/2 ; WHEN ABS(F) < EPS (SEE CONSTANTS, BELOW), SIN(F) = F ; OTHERWISE, THE ARGUMENT REDUCTION IS: ; LET X1 = [ABS(X)], THE GREATEST INTEGER IN ABS(X) ; X2 = ABS(X) - X1 ; N = THE NEAREST INTEGER TO ABS(X), FOR SIN, OR ; TO ABS(X) + PI/2, FOR COS ; THEN THE REDUCED ARGUMENT F = ((X1 - N*C1) + X2) - N*C2 ; WHERE C1 = 3.1104 (OCTAL) AND C2 IS GIVEN BELOW ; LET G = F**2 ; THEN R(G) = (G*XNUM/XDEN+RP1)*G ; WHERE XNUM = ((RP5*G+RP4)*G+RP3)*G+RP2 ; XDEN = ((G*Q2)*G+Q1)*G+Q0 ; AND RP5,RP4,RP3,RP2,RP1,Q2,Q1, AND Q0 ARE GIVEN BELOW ; AND SIN(F) = F + F*R(G). THE R(I) APPEAR BELOW ; FINALLY, ; SIN(X) = SIGN(X)*SIN(F)*(-1)**N, AND ; COS(X) = SIN(X + PI/2) ;THE RANGE OF DEFINITION FOR SIN IS ABS(X) <= 6746518852. ABS(X)+PI/2, IN ; COS(X), MUST BE LESS THAN 6746518852. SIN(X) = COS(X) = 0.0 AND AN ; ERROR MESSAGE WILL RESULT FOR ARGUMENTS OUT OF THAT RANGE ;REQUIRED (CALLED) ROUTINES: NONE ;REGISTERS T2, T3, T4, AND T5 WERE SAVED, USED, AND RESTORED ;THE ROUTINE HAS THE FOLLOWING CALLING SEQUENCE: ; XMOVEI L,ARG ; PUSHJ P,DSIN ; OR ; PUSHJ P,DCOS ;THE HIGH ORDER PART OF THE ANSWER IS RETURNED IN ACCUMULATOR T0. ;THE LOW ORDER PART OF THE ANSWER IS RETURNED IN ACCUMULATOR T1. SEARCH FORPRM TWOSEG 400000 SALL HELLO (DCOS,.) ;ENTRY TO COSINE ROUTINE DMOVE T0,@(L) ;OBTAIN ARGUMENT PUSH P,T2 ;SAVE AN ACCUMULATOR HRLZI T2,201400 ;SET SGN TO 1. JUMPGE T0,XPOS ;IF X IS NEGATIVE DMOVN T0,T0 ;Y = -X XPOS: DMOVEM T0,XDEN ;SAVE A COPY OF ABS(X) SETOM FLAG ;SET FLAG FOR COSINE ENTRY DFAD T0,PID2 ;Y = Y+PID2 JRST ALG ;PROCEED TO MAIN ALGORITHM HELLO (DSIN,.) ;ENTRY TO SIN ROUTINE DMOVE T0,@(L) ;GET ARG PUSH P,T2 ;SAVE AN ACCUMULATOR HRLZI T2,201400 ;SET SIGN TO 1.0 SETZM FLAG ;CLEAR FLAG FOR SINE ROUTINE JUMPGE T0,ALG ;IF ARG IS NEGATIVE DMOVN T0,T0 ;THEN, SET Y=-X AND HRLZI T2,576400 ;SET VALUE OF SIGN TO -1. ALG: MOVEM T2,SGN ;MOVE SIGN TO MEMORY CAMGE T0,YMAX1 ;IF HI OF Y HI OF YMAX JRST ERR1 ;GO TO ERR1 CAMGE T1,YMAX2 ;HI OF Y=HI OF YMAX, IS LO OF Y < LO OF YMAX? JRST YOK ;YES, PROCEED TO YOK ERR1: LERR (LIB,%,DSIN or DCOS: ABS(arg) too large; result = zero) MOVEI T0,0 ;SET RESULT MOVEI T1,0 ;TO ZERO POP P,T2 ;RESTORE ACCUMULATOR GOODBY (1) ; YOK: PUSH P,T3 ;SAVE ACCUMULATORS PUSH P,T4 ; PUSH P,T5 DMOVE T2,T0 ;SAVE A COPY OF ABS(X) SKIPE FLAG ;IF THIS IS THE COS ROUTINE DMOVE T2,XDEN ;RESTORE ABS(X) TO T2 DFMP T0,ODPI ;RN = Y/PI JFCL DMOVE T4,T0 ;SAVE A COPY OF RN CAMG T0,TWO26 ;IF RN IS .LE. 2**26 JRST LESS ;THEN GO TO LESS ASH T4,-33 ;SHIFT OFF MANTISSA SUBI T4,275 ;OBTAIN SHIFTING FACTOR ASH T1,(T4) ;SHIFT SECOND WORD APPROPRIATELY MOVE T5,T1 ;SAVE BITS OF SECOND WORD ASH T1,-1 ;SHIFT OFF REMAINDER MOVN T4,T4 ;PREPARE TO SHIFT BACK ADDI T4,1 ;ADD 1 TO SHIFTING FACTOR ASH T1,(T4) ;SHIFT BACK SECOND WORD TRNN T5,1 ;IF RIGHTMOST OF SAVED BITS IS OFF JRST NORND ;THEN GO TO NORND DFAD T0,ONE ;OTHERWISE ADD ONE TO XN ASH T5,-1 ;SHIFT REMAINDER OFF OF N ADDI T5,1 ;RND N UP JRST GOTN ;GO TO GOTN NORND: ASH T5,-1 ;SHIFT OFF FRACTIONAL PART JRST GOTN ;GO TO GOTN LESS: FIXR T0,T0 ;FIX RN MOVE T5,T0 ;SAVE N FLTR T0,T0 ;XN=FLOAT(N) MOVEI T1,0 ;ZERO SECOND WORD GOTN: TRNE T5,1 ;IS N ODD? MOVNS SGN ;YES,NEGATE SIGN SKIPE FLAG ;IF THE COSINE IS WANTED DFAD T0,PT5 ;THEN XN=XN-.5 DMOVE T4,T0 ;SAVE A COPY OF XN DFMP T0,C1 ;FORM XN*C1 DFSB T2,T0 ;GET |X| - XN * C1 DMOVE T0,T4 ;GET A COPY OF XN DFMP T4,C2 ;GET XN * C2 DFMP T0,C3 ;GET XN * C3 DFSB T2,T4 ;F = (|X| - XN*C1) - XN*C2 DFSB T2,T0 ; - XN*C3 DMOVE T0,T2 ;SAVE A COPY OF F JUMPGE T2,FPOS ;IF F IS NEGATIVE DMOVN T2,T2 ;F = -F FPOS: CAML T2,EPS ;IF F IS .GE. EPS JRST GEEPS ;GO TO GEEPS SKIPGE SGN ;IF SGN IS NEGATIVE DMOVN T0,T0 ;F = -F POP P,T5 ;RESTORE ACCUMULATORS POP P,T4 POP P,T3 POP P,T2 GOODBY (1) ;RETURN GEEPS: DFMP T2,T2 ;G = F*F DMOVE T4,T2 ;SAVE A COPY OF G DFAD T2,Q2 ;XDEN = G+Q2 DFMP T2,T4 ;*G DFAD T2,Q1 ;+Q1 DFMP T2,T4 ;*G DFAD T2,Q0 ;+Q0 DMOVEM T2,XDEN ;MOVE XDEN TO MEMORY DMOVE T2,T4 ;GET A COPY OF G DFMP T2,RP5 ;XNUM = G*RP5 DFAD T2,RP4 ;+RP4 DFMP T2,T4 ;*G DFAD T2,RP3 ;+RP3 DFMP T2,T4 ;*G DFAD T2,RP2 ;+RP2 DFMP T2,T4 ;*G DFDV T2,XDEN ;R(G) = XNUM/XDEN DFAD T2,RP1 ;+RP1 DFMP T2,T4 ;*G DFMP T2,T0 ;*F DFAD T0,T2 ;+F SKIPGE SGN ;IF SGN IS NEGATIVE DMOVN T0,T0 ;THEN NEGATE T0 RET: POP P,T5 ;RESTORE ACCUMULATORS POP P,T4 ; POP P,T3 POP P,T2 GOODBY (1) PID2: DOUBLE 201622077325,021026430215 ;PI/2 ODPI: DOUBLE 177505746033,162344202512 ;1/PI EPS: DOUBLE 142400000000,000000000000 ;2**(-31) C1: DOUBLE 202622077325,020000000000 ;HIGH 34 BITS OF PI C2: DOUBLE 140413214106,220000000000 ;NEXT 31 BITS OF PI C3: DOUBLE 101423063050,270033407151 ;C1+C2+C3=PI TO XTRA PREC Q2: DOUBLE 211612731352,270761316154 ;0.394924723520450141D+3 Q1: DOUBLE 221422322352,120271537275 ;0.702492288221842518D+5 Q0: DOUBLE 227512520255,264021153121 ;0.541748285645351853D+7 RP5: DOUBLE 605161526726,356444271113 ;-0.121560740596710190D-1 RP4: DOUBLE 203422023017,204642007402 ;0.428183075897778265D+01 RP3: DOUBLE 566026406450,010567611157 ;-0.489487151969463797D+03 RP2: DOUBLE 220540546606,025275050734 ;0.451456904704461990D+05 RP1: DOUBLE 601252525252,252525252421 ;-.166666666666666667D0 TWO26: DOUBLE 233400000000,000000000000 ;2**26 ONE: DOUBLE 201400000000,000000000000 ;1.0 PT5: DOUBLE 577400000000,000000000000 ;-.5 YMAX1: 241622077325 ;YMAX = (PI*2**31) YMAX2: 021026430215 RELOC ;DATA FLAG: 0 SGN: DOUBLE 0,0 XDEN: DOUBLE 0,0 RELOC PRGEND TITLE DSINH HYPERBOLIC SINE FUNCTION ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. JUNE 7, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DSINH EXTERN DSINH. DSINH=DSINH. PRGEND TITLE DSINH. HYPERBOLIC SINE FUNCTION ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. JUNE 7, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DSINH(X) IS CALCULATED AS FOLLOWS ; LET V BE APPROXIMATELY 2 SO THAT LNV AND ABS(X)+LN V ; CAN BE EXACTLY REPRESENTED WHEN X IS EXACTLY REPRESENTABLE. ; THEN, LETTING Y = ABS(X), AND NOTING THAT -SINH(-X)=SINH(X), ; FOR ; 0 <= Y < EPS, DSINH = Y*SIGN(X) ; EPS <= Y <= 1, DSINH = X-X*(Z*R(Z)), WHERE Z = Y**2 AND ; R(Z) IS GIVEN BELOW. ; 1 < Y <= 88.029678, DSINH = SIGN(X)*(EXP(Y)-EXP(-Y))/2, ; 88.029678 <= Y < 128 * LN(2) ; DSINH = SIGN(X)*((V/2)*EXP(Y - LN V)), ; Y >= 128 * LN(2), SINH = +MACHINE INFINITY * SIGN(X) ; LET Z = Y**2. THEN ; R(Z) = (RP0 + Z*(RP1 + Z*(RP2 + Z*RP3)))/(Q0 + Z*(Q1 + Z*(Q2 + Z))) ; WHERE RP0, RP1, RP2, Q0, Q1, AND Q2 ARE GIVEN BELOW. ;THE RANGE OF DEFINITION FOR DSINH IS ABS(X) < 128 * LN(2) AND ERROR MESSAGES ; WILL RESULT FOR ARGUMENTS OUT OF THAT RANGE. DSINH WILL BE SET ; TO + MACHINE INFINITY * SIGN(X). ;REQUIRED (CALLED) ROUTINES: DEXP ;REGISTERS T2, T3, T4, AND T5 WERE SAVED, USED, AND RESTORED ;THE ROUTINE HAS THE FOLLOWING CALLING SEQUENCE: ; XMOVEI L,ARG ; PUSHJ P,DSINH ;THE HIGH ORDER PART OF THE ANSWER IS RETURNED IN ACCUMULATOR T0 ;THE LOW ORDER PART OF THE ANSWER IS RETURNED IN ACCUMULATOR T1 SEARCH FORPRM TWOSEG 400000 SALL HELLO (DSINH,.) ;ENTRY TO DSINH ROUTINE DMOVE T0,@(L) ;OBTAIN X PUSH P,T5 ;SAVE REGISTER T5 HRLZI T5,201400 ;SET FLAG TO 1 JUMPGE T0,XPOS ;IF X IS NEGATIVE MOVN T5,T5 ;NEGATE FLAG DMOVN T0,T0 ;Y = -X XPOS: CAMLE T0,ONE ;IF Y IS .GT. 1.0 JRST DCSH ;GO TO DCSH LE: CAMG T0,TWOM31 ;IF Y IS .LE. 2**-31 JRST RET2 ;GO RETURN JUMPGE T5,NXT ;IF X IS NEGATIVE DMOVN T0,T0 ;Y = X NXT: PUSH P,T2 ;SAVE ACCUMULATORS PUSH P,T3 PUSH P,T4 DMOVE T2,T0 ;SAVE A COPY OF X DMOVNM T0,TEMP ;MOVE A COPY OF -X TO MEMORY DFMP T2,T2 ;Z = X*X DMOVE T4,T2 ;SAVE A COPY OF Z DFAD T4,Q2 ;XDEN = Z + Q2 DFMP T4,T2 ;*Z DFAD T4,Q1 ;+Q1 DFMP T4,T2 ;*Z DFAD T4,Q0 ;+Q0 DMOVEM T2,Z ;MOVE A COPY OF Z TO MEMORY DFMP T2,RP3 ;XNUM = Z*RP3 DFAD T2,RP2 ;+RP2 DFMP T2,Z ;*Z DFAD T2,RP1 ;+RP1 DFMP T2,Z ;*Z DFAD T2,RP0 ; +RP0 DFDV T2,T4 ;R(Z) = XNUM/XDEN DFMP T2,Z ;*Z DFMP T2,TEMP ;*(-X) DFAD T0,T2 ;+X POP P,T4 ;RESTORE ACCUMULATORS POP P,T3 POP P,T2 POP P,T5 GOODBY (1) ;RETURN DCSH: CAMGE T0,TWENT2 ;IF ABS(X) < 22, GO TO JRST ALG2 ; FULL CALCULATION CAMG T0,BIGX ;IF Y IS .LE. BIGX JRST EASY ;THEN GO TO EASY CAMGE T0,YMAXHI ;IF HI OF Y < YMAXHI JRST GETW ; GO TO GETW CAME T0,YMAXHI ;IF HI OF Y > YMAXHI JRST ERRD ;GO TO ERRD CAMGE T1,YMAXLO ;HI OF Y = YMAXHI JRST GETW ;IF LO OF Y .LE. YMAXLO, GO TO GETW ERRD: HRLOI T0,377777 ;SET RESULT TO SETO T1, ;MACHINE INFINITY JUMPG T5,DSNH ;IF X IS NEGATIVE DMOVN T0,T0 ;RESULT = -RESULT DSNH: LERR (LIB,%,DSINH: result overflow) POP P,T5 ;RESTORE ACCUMULATOR GOODBY (1) ;RETURN GETW: DFAD T0,LNV ;W = Y-LNV DMOVEM T0,TEMP ;MOVE W TO A TEMP REGISTER FUNCT DEXP., ;Z = EXP(W) DMOVEM T0,TEMP ;SAVE A COPY OF Z FMP T0,CON1 ;RESULT = Z*CON1 SETZ T1, ;ZERO SECOND WORD DFAD T0,TEMP ;+Z JFCL ERRD ;OVERFLOW POSSIBLE RET2: JUMPGE T5,RET3 ;IF SGNFLG IS NEGATIVE DMOVN T0,T0 ;RESULT = -RESULT RET3: POP P,T5 ;RESTORE ACCUMULATOR GOODBY (1) ;RETURN EASY: DMOVEM T0,TEMP ;MOVE Y TO TEMP FUNCT DEXP., ;GET ITS EXP FSC T0,-1 ;DIV BY 2 JRST RET2 ;GET RIGHT SIGN FOR RESULT. ALG2: PUSH P,T2 ;SAVE MORE ACCUMULATORS PUSH P,T3 PUSH P,T4 DMOVEM T0,TEMP ;MOVE Y TO A TEMPORARY REGISTER FUNCT DEXP., ;Z = DEXP(Y) DMOVN T2,T0 ;SAVE A COPY OF Z MOVEM T5,SGNFLG ;MOVE FLAG TO MEMORY HRLZI T4,201400 ;SET T4 TO 1.0 SETZ T5, ;CLEAR SECOND WORD DFDV T4,T2 ;1/Z DFAD T0,T4 ;+Z FSC T0,-1 ;*1/2 POP P,T4 ;RESTORE ACCUMULATORS POP P,T3 POP P,T2 RET1: SKIPGE SGNFLG ;IF SGNFLG IS NEGATIVE DMOVN T0,T0 ;RESULT = -RESULT RET: POP P,T5 ;RESTORE ACCUMULATOR GOODBY (1) ;RETURN RP0: DOUBLE 223527442325,224632030652 ;.35181283430177117881D+6 RP1: DOUBLE 216551270255,245012217565 ;.11563521196851768270D+5 RP2: DOUBLE 210507410130,341337112321 ;.16375798202630751372D+3 RP3: DOUBLE 200624234756,033744032643 ;.78966127417357099479D+0 Q0: DOUBLE 551376246137,720314355210 ;-.21108770058106271242D+7 Q1: DOUBLE 220432412710,355457306744 ;.36162723109421836460D+5 Q2: DOUBLE 566352207437,615500015770 ;-.27773523119650701667D+3 ONE: 201400000000 ;1.0 BIGX: 207540074635 ;88.029691 YMAXHI: 207542710277 ;88.722839111672999604 YMAXLO: 276434757153 ; TWOM31: 142400000000 ;2**-31 LNV: DOUBLE 577235067500,101343000000 ;-.693147180559947174D0 CON1: 120414520000 ;.186417721E-14 TWENT2: 205540000000 ;22 RELOC ;DATA SGNFLG: 0 TEMP: DOUBLE 0,0 Z: DOUBLE 0,0 RELOC PRGEND TITLE DSQRT SQUARE ROOT FUNCTION ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. MARCH 19, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DSQRT EXTERN DSQRT. DSQRT=DSQRT. PRGEND TITLE DSQRT. SQUARE ROOT FUNCTION ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. MARCH 19, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DSQRT(X) IS CALCULATED AS FOLLOWS ;THE ROUTINE CALCULATES THE SQUARE ROOT OF THE ABSOLUTE VALUE OF A ;DOUBLE PRECISION ARGUMENT BY DOING A LINEAR SINGLE PRECISION ;APPROXIMATION ON THE HIGH ORDER WORD, FOLLOWED BY TWO SINGLE PRECISION ;NEWTON ITERATIONS AND TWO DOUBLE PRECISION NEWTON ITERATIONS. AN ERROR ;MESSAGE RESULTS FOR NEGATIVE ARGUMENTS. ;REQUIRED (CALLED) ROUTINES: NONE ;REGISTERS T2, T3, T4, AND T5 WERE SAVED, USED, AND RESTORED ;THE ROUTINE HAS THE FOLLOWING CALLING SEQUENCE: ; XMOVEI L,ARG ; PUSHJ P,DSQRT ;THE HIGH ORDER PART OF THE ANSWER IS RETURNED IN ACCUMULATOR T0 ;THE LOW ORDER PART OF THE ANSWER IS RETURNED IN ACCUMULATOR T1 SEARCH FORPRM TWOSEG 400000 SALL HELLO (DSQRT,.) ;ENTRY TO DSQRT ROUTINE DMOVE T0,@(L) ;GET DP ARGUMENT JUMPG T0,DSQRTP ;ARGUMENT IS GREATER THAN 0 JUMPE T0,DSQRT4 ;ARGUMENT IS ZERO LERR (LIB,%,DSQRT: negative arg; result = SQRT(ABS(arg))) DMOVN T0,T0 ;ARG = -ARG DSQRTP: PUSH P,T2 ;SAVE ACCUMULATORS PUSH P,T3 PUSH P,T4 PUSH P,T5 MOVE T2,T0 ;GET SPARE COPY OF HIGH ORDER BITS MOVE T3,T0 ;GET SPARE COPY OF HIGH ORDER BITS MOVE T4,T0 ;GET SPARE COPY OF HIGH ORDER BITS LSH T3,-1 TLZE T3,400 ;WAS EXPONENT ODD? JRST DSQRT2 ;YES, GO TO DSQRT2 ADD T3,[XWD 267,607000] ;NO,COMPUTE LINEAR APPROXIMATION FMPRI T3,301454 ;RESCALE JRST DSQRT3 DSQRT2: ADD T3,[XWD 267,607000] ;ODD EXPONENT, COMPUTE LINEAR APPROXIMATION FMPRI T3,301650 ;RESCALE DSQRT3: FDV T4,T3 ;ORIGINAL ARG / INITIAL GUESS FAD T4,T3 ;AVERAGE THEM FSC T4,-1 ; FDV T2,T4 ;ORIGINAL ARG / SECOND GUESS FAD T2,T4 ;AVERAGE THEM FSC T2,-1 MOVEI T3,0 ;LOW HALF OF 1ST APPROX. IS 0 DMOVE T4,T0 ;GET ARG BACK DFDV T4,T2 ;ORIGINAL ARG / THIRD GUESS DFAD T2,T4 ;AVERAGE THEM FSC T2,-1 ; DFDV T0,T2 ;ORIGINAL ARG / FOURTH GUESS DFAD T0,T2 ;AVERAGE THEM FSC T0,-1 POP P,T5 ;RESTORE ACCUMULATORS POP P,T4 POP P,T3 POP P,T2 DSQRT4: GOODBY (1) ;RETURN PRGEND TITLE DCOTAN TANGENT AND COTANGENT FUNCTIONS ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. APRIL 16, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DCOTAN EXTERN DCOTN. DCOTAN=DCOTN. PRGEND TITLE DTAN TANGENT AND COTANGENT FUNCTIONS ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC APRIL 16, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DTAN EXTERN DTAN. DTAN=DTAN. PRGEND TITLE DTAN. TANGENT AND COTANGENT FUNCTIONS ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC APRIL 16, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;TAN(X) AND COTAN(X) ARE CALCULATED AS FOLLOWS; ; THE IDENTITIES ; TAN(PI/2.0-G)=1.0/TAN(G) ; TAN(N*PI+H)=TAN(H), -PI/2.0 < H <= PI/2.0 ; TAN(-X)=-TAN(X) ; COTAN(X)=1.0/TAN(X) ; COTAN(-X)=-COTAN(X) ; ARE USED TO REDUCE THE ORIGINAL PROBLEM TO A PROBLEM WITH ; AN ARGUMENT GREATER THAN -PI/2.0 AND LESS THAN OR EQUAL TO ; PI/2.0, PI=3.14159265358979324. ; THEN DEFINE N AND F SO THAT ; X=N*PI/4.0+F, 0.0 <= F <= PI/4.0, WITH PI=3.14159265358979324. ; THEN ; TAN(F) = F, IF F<2**(-31) ; = R(F), OTHERWISE ; WHERE ; R(F)=(((XP3*G+XP2)*G+XP1)*G+XP0)/((((Q4*G+Q3)*G+Q2)*G+Q1)*G+Q0) ; AND ; G = F*F ; XPi, Qi ARE GIVEN BELOW ; THE RESULT IS THEN RECIPROCATED, IF NECESSARY, AND GIVEN ; THE APPROPRIATE SIGN. ; THE APPROXIMATION USED IS A MODIFICATION OF HART 4287. ;THE RANGE OF DEFINITION FOR DTAN IS 0 < ABS(X) < ((2**31) * (PI/2)) = 3373259426.D0 ; AND FOR DCOTAN(X), (2**(-127))*(1+(2**(-61))) < ABS(X) < ((2**31)*(PI/2)) IS NECESSARY. ; OTHERWISE, ERROR MESSAGES WILL RESULT. ;REQUIRED (CALLED) ROUTINES: NONE ;REGISTERS T2, T3, T4, AND T5 WERE SAVED, USED, AND RESTORED ;THE ROUTINE HAS THE FOLLOWING CALLING SEQUENCE: ; XMOVEI L,ARG ; PUSHJ P,DTAN ; OR ; PUSHJ P,DCOTAN ;THE HIGH ORDER PART OF THE ANSWER IS RETURNED IN ACCUMULATOR T0, ;THE LOW ORDER PART OF THE ANSWER IS RETURNED IN ACCUMULATOR T1. SEARCH FORPRM TWOSEG 400000 SALL HELLO (DCOTN.,DCOTAN) ;ENTRY TO DCOTAN ROUTINE PUSH P,T4 ;SAVE ACCUMULATORS PUSH P,T5 SETZM FLAG ;CLEAR FLAG FOR DCOTAN DMOVE T0,@(L) ;GET ARG DMOVE T4,T0 ;Y = X JUMPGE T0,XPOS ;IF X IS NEGATIVE DMOVN T0,T0 ;SET Y = -X XPOS: CAMGE T0,EPS1 ;IF THE HI ORDER PART OF Y IS TOO SMALL JRST ERR1 ;THEN GO TO ERR1 CAME T0,EPS1 ;IF HI OF Y > HI OF EPS JRST YOK ;THEN GO TO YOK JUMPN T1,YOK ;HI OF Y = HI OF EPS, IS LO OF Y OK? ERR1: LERR (LIB,%,DCOTAN: ABS(arg) too large; result = zero) HRLOI T0,377777 ;SET RESULT TO MACHINE HRLOI T1,377777 ;INFINITY JUMPGE T4,RET ;IF ARG IS NEGATIVE DMOVN T0,T0 ;RESULT = - RESULT POP P,T5 ;RESTORE ACCUMULATORS POP P,T4 GOODBY (1) ;RETURN HELLO (DTAN,.) ;ENTRY TO TAN ROUTINE PUSH P,T4 ;SAVE ACCUMULATORS PUSH P,T5 SETOM FLAG ;SET FLAG FOR DTAN ENTRY DMOVE T0,@(L) ;GET COPY OF ARG DMOVE T4,T0 ;Y = X JUMPGE T0,YOK ;IF X IS NEGATIVE DMOVN T0,T0 ;Y = -X YOK: CAMGE T0,YMAX1 ;IF HI OF Y < HI OF YMAX JRST ALG ;PROCEED TO ALG CAME T0,YMAX1 ;IF HI OF Y > HI OF YMAX JRST ERR2 ;GO TO ERR2 CAMG T1,YMAX2 ;HI OF Y=HI OF YMAX, IS LO OF Y .LE. LO OF YMAX? JRST ALG ;YES, PROCEED TO ALG ERR2: LERR (LIB,%,DTAN or DCOTAN: result underflow) MOVEI T0,0 ;SET RESULT TO ZERO MOVEI T1,0 ;SET SECOND WORD TO ZERO POP P,T5 ;RESTORE ACCUMULATORS POP P,T4 GOODBY (1) ;RETURN ALG: PUSH P,T2 ;SAVE ACCUMULATORS PUSH P,T3 CAML T0,PIO4HI ;COMPARE |ARG| WITH PI/4 JRST REDUCE ;REDUCTION IS NECESSARY DMOVE T2,T0 ;|F| TO T0 DMOVE T0,T4 ;AND SIGNED F TO T0. SETZM N ;STORE 0 FOR N. JRST FPOS ;BYPASS REDUCTION REDUCE: DFMP T0,TWODPI ;Y*(2/PI) MOVE T2,T0 ;SAVE A COPY OF Y*(2/PI) CAMG T2,TWO26 ;IF Y*(2/PI) IS .LE. 2**26 JRST LES ;THEN GO TO LES ASH T2,-33 ;SHIFT OFF MANTISSA SUBI T2,275 ; ASH T1,(T2) ;SHIFT SECOND WORD APPROPRIATELY MOVE T3,T1 ;SAVE BITS OF SECOND WORD ASH T1,-1 ;SHIFT REMAINDER OFF MOVN T2,T2 ;PREPARE TO SHIFT BACK ADDI T2,1 ; ASH T1,(T2) ;SHIFT BACK TRNN T3,1 ;IF RIGHTMOST OF SAVED BITS IS OFF JRST NORND ;THEN GO TO NORND DFAD T0,ONE ;OTHERWISE, ADD 1 TO XN ASH T3,-1 ;PREPARE TO LOOK AT NEXT BIT ADDI T3,1 ;ADD 1 TO N JRST RND ;GO TO RND NORND: ASH T3,-1 ;SHIFT OFF FRACTIONAL PART RND: MOVEM T3,N ;MOVE N TO MEMORY JRST GOTIT ;GO TO GOTIT LES: FIXR T0,T0 ;FIX IT MOVEM T0,N ;MOVE N TO MEMORY FLTR T0,T0 ;FLOAT IT MOVEI T1,0 ;ZERO SECOND WORD GOTIT: JUMPLE T4,NEXT ;IF X IS POSITIVE DMOVN T0,T0 ;NEGATE XN NEXT: DMOVE T2,T0 ;GET A COPY OF -XN DFMP T2,C1 ;FORM -XN*C1 DFAD T4,T2 ;GET |X| - N*C1 DMOVE T2,T0 ;GET A COPY OF -XN DFMP T0,C2 ;GET -XN*C2 DFMP T2,C3 ;GET -XN*C3 DFAD T2,T4 ;F = (|X| - XN*C1) - XN*C2 DFAD T2,T0 ; - XN*C3. DMOVE T0,T2 ;MOVE COPY OF F INTO T0 JUMPGE T2,FPOS ;IF F IS NEGATIVE DMOVN T2,T2 ;F = -F FPOS: CAML T2,HIEPS ;IS F < EPS? JRST NO ;NO, GO TO NO HRLZI T4,201400 ;YES , F< EPS, MOVEI T5,0 ;XDEN = 1.0 JRST FLGCHK ;GO TO CHECK FLAG NO: DFMP T2,T2 ;NO, F .GE. EPS, SET G=F*F DMOVE T4,T2 ;SAVE A COPY OF G DFMP T2,XP4 ;XNUM = XP4 * G + DFAD T2,XP3 ; P3 * DFMP T2,T4 ; G + DFAD T2,XP2 ; P2 * DFMP T2,T4 ; G + DFAD T2,XP1 ; P1 * DFMP T2,T4 ; G * DFMP T2,T0 ; F + DFAD T0,T2 ; F DMOVE T2,T4 ;SAVE A COPY OF G DFMP T4,Q4 ;XDEN = Q4 * G + DFAD T4,Q3 ; Q3 * DFMP T4,T2 ; G + DFAD T4,Q2 ; Q2 * DFMP T4,T2 ; G + DFAD T4,Q1 ; Q1 * DFMP T4,T2 ; G + DFAD T4,ONE ; 1.0 FLGCHK: MOVE T2,N ;GET COPY OF N SKIPE FLAG ;IF FLAG IS NOT ZERO JRST DA ;THEN GO TO DA TRNN T2,1 ;OTHERWISE, IS N ODD? JRST DOVN ;NO, GO GET RESULT DMOVN T0,T0 ;XNUM = -XNUM JRST NOVD ;GO TO NOVD DA: TRNN T2,1 ;FLAG MUST BE ONE, IS N ODD? JRST NOVD ;NO, GO TO NOVD DMOVN T0,T0 ;XNUM = -XNUM DOVN: DFDV T4,T0 ;RESULT = XDEN/XNUM DMOVE T0,T4 JRST RET ;GO TO RET NOVD: DFDV T0,T4 ;RESULT = XNUM/XDEN RET: POP P,T3 ;RESTORE ACCUMULATORS POP P,T2 POP P,T5 POP P,T4 GOODBY (1) XP1: DOUBLE 601346652066,115432245623 ;-.1372889460941120802 XP2: DOUBLE 171401224404,126140152543 ; .3925934686364577602E-02 XP3: DOUBLE 616034314474,026356056117 ;-.2882482747560198194E-04 XP4: DOUBLE 147766720604,360442362056 ; .2927308283322907641E-07 Q1: DOUBLE 600036052305,321342375437 ;-.4706222794274454135 Q2: DOUBLE 173702007252,226306364361 ; .2746669449551304872E-01 Q3: DOUBLE 612131325464,137044675153 ;-.4030063705745304384E-03 Q4: DOUBLE 155540343710,045637644377 ; .1312960309685759549E-05 C1: DOUBLE 201622077325,020000000000 ;FIRST 34 BITS OF PI/2 C2: DOUBLE 137413214106,220000000000 ;NEXT 31 BITS OF PI/2 C3: DOUBLE 100423063050,270033407151 ;NEXT 62 BITS OF PI/2. EPS1: 002400000000 ;HIGH ORDER PART OF EPS1 YMAX1: 240622077325 ;HIGH ORDER PART OF YMAX YMAX2: 021026430215 ;LOW ORDER PART OF YMAX HIEPS: 142400000000 ; 2**(-31) PIO4HI: 200622077325 ;HIGH ORDER PART OF PI/4 TWODPI: DOUBLE 200505746033,162344202512 ONE: DOUBLE 201400000000,000000000000 TWO26: 233400000000 ;2**26 RELOC ;DATA N: 0 FLAG: 0 RELOC PRGEND TITLE DTANH HYPERBOLIC TANGENT FUNCTION ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. JUNE 1, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DTANH EXTERN DTANH. DTANH=DTANH. PRGEND TITLE DTANH. HYPERBOLIC TANGENT FUNCTION ; (DOUBLE PRECISION FLOATING POINT) SUBTTL IMSL, INC. JUNE 1, 1979 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. ;TANH(X) IS CALCULATED AS FOLLOWS ; TANH(X) = -TANH(-X). LET F = ABS(X) ; IF F > 22.1807100, TANH = 1.0*SIGN(X) ; IF F > LN(3)/2 AND F <= 22.1807100, TANH = RESULT 1 = ; SIGN(X)*(1 - 2/(EXP(2*F) + 1))) ; IF F < 2**(-31)*SQRT(3), TANH = F*SIGN(X) ; OTHERWISE, LET G = F**2, AND OBTAIN RESULT 2 AS ; TANH = SIGN(X)*(F + F*R(G)), WHERE ; R(G) = G*(RP0+G*(RP1+RP2*G))/(Q0+G*(Q1+G*(Q2+G))). ; RP0, RP1, RP2, Q0, Q1, AND Q2 APPEAR BELOW. ;THE RANGE OF DEFINITION FOR TANH IS THE REPRESENTABLE REAL NUMBERS ;REQUIRED (CALLED) ROUTINES: DEXP ;REGISTERS T2, T3, T4, AND T5 WERE SAVED, USED, AND RESTORED ;THE ROUTINE HAS THE FOLLOWING CALLING SEQUENCE: ; XMOVEI L,ARG ; PUSHJ P,TANH ;THE HIGH ORDER PART OF THE ANSWER IS RETURNED IN ACCUMULATOR T0 ;THE LOW ORDER PART OF THE ANSWER IS RETURNED IN ACCUMULATOR T1 SEARCH FORPRM TWOSEG 400000 SALL HELLO (DTANH,.) ;ENTRY TO DTANH ROUTINE DMOVE T0,@(L) ;OBTAIN X PUSH P,T5 ;SAVE AN ACCUMULATOR MOVE T5,T0 ;SAVE A COPY OF X JUMPGE T0,XPOS ;IF X IS NEGATIVE DMOVN T0,T0 ;F = -X XPOS: CAMG T0,LN2T32 ;IF F IS .LE. LN2T32 JRST CALCF ;GO TO CALCF HRLZI T0,201400 ;SET RESULT TO 1.0 SETZ T1, ;ZERO SECOND WORD JUMPGE T5,RET1 ;IF X IS NEGATIVE DMOVN T0,T0 ;RESULT = -RESULT RET1: POP P,T5 ;RESTORE ACCUMULATORS GOODBY (1) ;RETURN CALCF: PUSH P,T2 ;SAVE MORE ACCUMULATORS PUSH P,T3 PUSH P,T4 CAMG T0,LN3D2 ;IF F IS .LE. LN3D2 JRST ALG1 ;GO TO ALG1 FSC T0,1 ;F = F+F DMOVEM T0,TEMP ;SAVE F IN A TEMP REGISTER FUNCT DEXP., ;EXP(2*F) DFAD T0,ONE ;+1.0 HRLZI T2,575400 ;SET T2 TO -2.0 SETZ T3, ;ZERO SECOND WORD DFDV T2,T0 ;-2/(EXP(2*F)+1) DMOVE T0,T2 DFAD T0,ONE ;1 - 2/(EXP(2*F)+1) JUMPGE T5,RET ;IF X IS NEGATIVE DMOVN T0,T0 ;RESULT = -RESULT POP P,T4 ;RESTORE ACCUMULATORS POP P,T3 POP P,T2 POP P,T5 GOODBY (1) ;RETURN ALG1: CAML T0,CON1 ;IF F IS .GE. .806549008E-9 JRST ALG2 ;GO TO ALG2 JUMPGE T5,RET ;IF X IS NEGATIVE DMOVN T0,T0 ;RESULT = -RESULT POP P,T4 ;RESTORE ACCUMULATORS POP P,T3 POP P,T2 POP P,T5 GOODBY (1) ;RETURN ALG2: JUMPGE T5,POSARG ;IF ARGUMENT < 0 DMOVN T0,T0 ; NEGATE |ARGUMENT| POSARG: DMOVEM T0,TEMP ;SAVE SIGNED ARGUMENT DFMP T0,T0 ;G = F*F DMOVE T2,T0 ;SAVE A COPY OF G DFAD T2,Q2 ;XDEN = G+Q2 DFMP T2,T0 ;*G DFAD T2,Q1 ;+Q1 DFMP T2,T0 ;*G DFAD T2,Q0 ; + Q0 DMOVE T4,T0 ;SAVE A COPY OF G DFMP T4,RP2 ;XNUM = G*RP2 DFAD T4,RP1 ; + RP1 DFMP T4,T0 ; * G DFAD T4,RP0 ; + RP0 DFMP T0,T4 ; *G DFDV T0,T2 ;R(G) = XNUM/XDEN DFMP T0,TEMP ;RESULT = F*R DFAD T0,TEMP ; + F RET: POP P,T4 ;RESTORE ACCUMULATORS POP P,T3 POP P,T2 POP P,T5 GOODBY (1) ;RETURN RP0: DOUBLE 564154513215,220360747434 ;-.161341190239962281D+4 RP1: DOUBLE 570163061227,221321463447 ;-.992259296722360833D+2 RP2: DOUBLE 577022172714,101054201376 ;-.964374927772254698D0 Q0: DOUBLE 215456407425,323513222652 ;.484023570719886887D+4 Q1: DOUBLE 214427161323,100370362574 ;.22337720718962312926D+4 Q2: DOUBLE 207702765170,172773772374 ;.112744743805349493D+3 LN2T32: 205542710300 ;22.1807100E0 CON1: 141673317272 ;2**-32 * SQRT(3) LN3D2: 200431175236 ;.549306144E0 HALF: DOUBLE 200400000000,000000000000 ;1/2 ONE: DOUBLE 201400000000,000000000000 ;1.0 RELOC ;DATA TEMP: DOUBLE 0,0 RELOC PRGEND COMMENT | Not in version 6 TITLE DINT DOUBLE PRECISION TRUNCATION TO INTEGER SUBTTL CHRIS SMITH/CKS 29-Jan-80 ;COPYRIGHT (C) 1980,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DINT EXTERN DINT. DINT=DINT. PRGEND TITLE DINT. DOUBLE PRECISION TRUNCATION TO INTEGER SUBTTL CHRIS SMITH/CKS 29-Jan-80 ;COPYRIGHT (C) 1980,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM TWOSEG 400000 SALL HELLO (DINT,.) PUSH P,T2 ;SAVE TEMP REGISTERS PUSH P,T3 PUSH P,T4 DMOVE T0,@(L) ;GET ABS(ARG) IN T0-T1 JUMPGE T0,DINT1 DMOVN T0,T0 DINT1: CAMGE T0,[201400000000] ;IS NUMBER .LT. 1? JRST [SETZB T0,T1 ;YES, DINT IS 0 JRST R1] CAML T0,[277400000000] ;IS NUMBER ALL INTEGER? JRST R ;YES, FAST RETURN LDB T4,[POINT 8,T0,8] ;GET EXPONENT TRC T4,-1 ;GET ONE'S COMPLEMENT FOR RIGHT SHIFT MOVSI T2,777000 ;GET INITIAL MASK TO GET RID OF FRACTION SETZ T3, ASHC T2,201(T4) ;MOVE MASK OVER INTEGER PART AND T0,T2 ;CLEAR FRACTION BITS AND T1,T3 R: SKIPGE @(L) ;ATTACH ORIGINAL SIGN DMOVN T0,T0 R1: POP P,T4 ;RESTORE TEMP REGISTERS POP P,T3 POP P,T2 POPJ P, ;RETURN PRGEND End of code not in version 6 | TITLE DABS DOUBLE PRECISION ABSOLUTE VALUE SUBTTL H. P. WEISS/HPW 11-DEC-73 ;COPYRIGHT (C) 1973,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DABS EXTERN DABS. DABS=DABS. PRGEND TITLE DABS. DOUBLE PRECISION ABSOLUTE VALUE SUBTTL /KK/DMN/TWE ;COPYRIGHT (C) 1980,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. ;DOUBLE PRECISION ABSOLUTE VALUE FUNCTION ;THIS ROUTINE RETURNS THE ABSOLUTE VALUE OF A DOUBLE PRECISION ;ARGUMENT SEARCH FORPRM TWOSEG 400000 SALL HELLO (DABS,.) ;[235] ENTRY TO DABS ROUTINE DMOVE T0,@(L) ;GET ARGUMENT JUMPGE T0,DABSX ;IF POSITIVE, RETURN DMOVN T0,T0 ;ELSE NEGATE IT DABSX: GOODBY (1) ;RETURN PRGEND TITLE DMAX1 DOUBLE PRECISION MAXIMUM OF A SERIES OF ARGUMENTS SUBTTL H. P. WEISS/HPW 11-DEC-73 ;COPYRIGHT (C) 1973,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DMAX1 EXTERN DMAX1. DMAX1=DMAX1. PRGEND TITLE DMAX1. DOUBLE PRECISION MAXIMUM OF A SERIES OF ARGUMENTS SUBTTL /TWE/CKS 29-Jan-80 ;COPYRIGHT (C) 1980,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM TWOSEG 400000 SALL HELLO (DMAX1,.) ;[235] ENTRY TO DMAX1 ROUTINE PUSH P,T2 ;SAVE AC T2 PUSH P,T3 ;Save AC t3 HLLZ T3,-1(L) ;Get the F10 arg count DMOVE T0,@(L) ;GET FIRST ARGUMENT JRST DMAX.2 ;ADDRESS OF NEXT, START CHECKING DMAX.1: XMOVEI T2,@0(L) ;GET ADDRESS OF NEXT ARGUMENT CAMLE T0,(T2) ;IS HIGH ORDER .GT. THIS ONE? JRST DMAX.2 ;YES, GET ADDRESS OF NEXT IN LIST CAME T0,(T2) ;ARE HIGH ORDER WORDS EQUAL? JRST DMAX.3 ;NO, REPLACEMENT NECESSARY CAMGE T1,1(T2) ;SKIP IF PRESENT ARG IS LARGER DMAX.3: DMOVE T0,(T2) ;PICK UP NEW ARG DMAX.2: ADDI L,1 ;Point to next arg AOBJN T3,DMAX.1 ;Loop for all args. POP P,T3 ;End of arg list POP P,T2 ;, restore acs GOODBY (2) ;Return answer in 0 and 1 PRGEND TITLE DMIN1 DOUBLE PRECISION MINIMUM OF A SERIES OF ARGUMENTS SUBTTL H. P. WEISS/HPW 11-DEC-73 ;COPYRIGHT (C) 1973,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DMIN1 EXTERN DMIN1. DMIN1=DMIN1. PRGEND TITLE DMIN1. DOUBLE PRECISION MINIMUM OF A SERIES OF ARGUMENTS SUBTTL /TWE/CKS 29-Jan-80 ;COPYRIGHT (C) 1980,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM TWOSEG 400000 SALL HELLO (DMIN1,.) ;[235] ENTRY TO DMIN1 ROUTINE PUSH P,T2 ;SAVE AC T2 PUSH P,T3 ;And T3 HLLZ T3,-1(L) ;Get the F10 arg count DMOVE T0,@(L) ;GET FIRST ARGUMENT JRST DMIN.2 ;ADDRESS OF NEXT, START CHECKING DMIN.1: XMOVEI T2,@0(L) ;GET ADDRESS OF NEXT ARGUMENT CAMGE T0,(T2) ;IS HIGH ORDER .LT. THIS ONE? JRST DMIN.2 ;YES, GET ADDRESS OF NEXT IN LIST CAME T0,(T2) ;ARE HIGH ORDER WORDS EQUAL? JRST DMIN.3 ;NO, REPLACEMENT NECESSARY CAMLE T1,1(T2) ;SKIP IF PRESENT ARG IS LARGER DMIN.3: DMOVE T0,(T2) ;PICK UP NEW ARG DMIN.2: ADDI L,1 ;Point to next arg AOBJN T3,DMIN.1 ;End of arg list? POP P,T3 ;Yes, restore T3 POP P,T2 ; and T2 GOODBY (2) ;Return answer in 0 and 1 PRGEND TITLE DSIGN DOUBLE PRECISION SIGN ROUTINE SUBTTL H. P. WEISS/HPW 11-DEC-73 ;COPYRIGHT (C) 1973,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DSIGN EXTERN DSIGN. DSIGN=DSIGN. PRGEND TITLE DSIGN. DOUBLE PRECISION SIGN ROUTINE SUBTTL /KK/DMN/CKS 29-Jan-80 ;COPYRIGHT (C) 1980,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. ;FROM LIB40 V.006. ;DOUBLE PRECISION TRANSFER OF SIGN ;THIS ROUTINE RETURNS ABSF(ARG1)*SIGN(ARG2) ;THE ROUTINE MAKES USE OF THE FOLLOWING TABLE: ;ARG1 ARG2 RESULT CHANGE OF SIGN? ;+ + + NO ;+ - - YES ;- + + YES ;- - - NO SEARCH FORPRM TWOSEG 400000 SALL HELLO (DSIGN,.) ;[235] ENTRY TO DSIGN ROUTINE DMOVE T0,@(L) ;GET FIRST ARGUMENT SKIPGE @1(L) ;THEN JUMPL T0,OUT ;CHOOSE SKIPL @1(L) ;THE JUMPGE T0,OUT ;CORRECT DMOVN T0,T0 ;SIGN. OUT: GOODBY (2) ;EXIT PRGEND TITLE DFLOAT INTEGER TO DOUBLE PRECISION CONVERSION SUBTTL H. P. WEISS/HPW 11-DEC-73 ;COPYRIGHT (C) 1973,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DFLOAT EXTERN DFLOT. DFLOAT=DFLOT. PRGEND TITLE DFLOT. INTEGER TO DOUBLE PRECISION CONVERSION SUBTTL ED YOURDON/KK/VAA/TWE/DMN ;COPYRIGHT (C) 1973,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM EXTERN DFL.0 NOSYM TWOSEG 400000 SALL HELLO (DFLOT.,DFLOAT) ;[235] ENTRY TO DFLOAT ROUTINE MOVE T0,@(L) ;PICK UP THE ARGUMENT PJRST DFL.0 ;USE DFL.0 ROUTINE PRGEND TITLE IDINT DOUBLE PRECISION TO INTEGER CONVERSION SUBTTL H. P. WEISS/HPW 11-DEC-73 ;COPYRIGHT (C) 1973,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY IDINT EXTERN IDINT. IDINT=IDINT. PRGEND TITLE IDINT. DOUBLE PRECISION TO INTEGER CONVERSION SUBTTL TOM EGGERS/DMN ;COPYRIGHT (C) 1973,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM EXTERN IDF.0 NOSYM TWOSEG 400000 SALL HELLO (IDINT,.) ;[235] ENTRY TO IDINT DFIX1: DMOVE T0,@(L) ;GET THE ARGUMENT PJRST IDF.0 ;USE IDF.0 ROUTINE PRGEND TITLE DBLE REAL TO DOUBLE PRECISION CONVERSION SUBTTL H. P. WEISS/HPW 11-DEC-73 ;COPYRIGHT (C) 1973,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY DBLE EXTERN DBLE. DBLE=DBLE. PRGEND TITLE DBLE. REAL TO DOUBLE PRECISION CONVERSION SUBTTL ED YOURDON/KK ;COPYRIGHT (C) 1973,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. ;FROM LIB40 V.005. ;SINGLE PRECISION TO DOUBLE PRECISION CONVERSION ;THIS ROUTINE CONVERTS A SINGLE PRECISION ARGUMENT TO A ;DOUBLE PRECISION ANSWER WITH THE LOW ORDER WORD SET TO 0 SEARCH FORPRM NOSYM TWOSEG 400000 SALL HELLO (DBLE,.) ;[235] ENTRY TO DOUBLE ROUTINE MOVE T0,@(L) ;GET THE ARGUMENT IN HIGH ORDER SETZ T1, ;SET LOW ORDER PORTION TO ZERO GOODBY (1) ;EXIT PRGEND TITLE SNGL DOUBLE PRECISION TO REAL CONVERSION SUBTTL H. P. WEISS/HPW 11-DEC-73 ;COPYRIGHT (C) 1973,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM NOSYM ENTRY SNGL EXTERN SNGL. SNGL=SNGL. PRGEND TITLE SNGL. DOUBLE PRECISION TO REAL CONVERSION SUBTTL /DMN/TWE ;COPYRIGHT (C) 1973,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. ;DOUBLE PRECISION TO SINGLE PRECISION CONVERSION FUNCTION ;THIS ROUTINE OBTAINS THE MOST SIGNIFICANT PART OF A DOUBLE ;PRECISION ARGUMENT. SEARCH FORPRM EXTERN SNG.0 NOSYM TWOSEG 400000 SALL HELLO (SNGL,.) ;[235] ENTRY TO SNGL ROUTINE DMOVE T0,@(L) ;GET THE ARGUMENT IN AC 0 PJRST SNG.0 ;USE AC0 SINGLE ROUTINE PRGEND UNIVERSAL FORDAR GENERAL DOUBLE PRECISION ARITHMETICS SUBTTL TOM EGGERS/DMN/DRT/JNG/SWG 14-Aug-79 ;COPYRIGHT (C) 1972,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SUBTTL MACROS FOR FORDAR UNIVERSAL FILE IF1,< ;ONLY ONCE ; DOUBLE PRECISION FLOAT FUNCTION "DFLOAT" DEFINE DFL (X) < XALL ENTRY DFL.'X ;ENTRY POINT TO DFL.'X SIXBIT /DFL.'X/ DFL.'X: MOVEI X+1,0 ;CLEAR LOW ORDER WORD ASHC X,-8 ;MAKE ROOM FOR EXPONENT IN HI WORD TLC X,243000 ;SET EXP TO 27+8 DECIMAL DFAD X,[EXP 0,0] ;NORMALIZE POPJ P, ;RETURN X=THE DOUBLE PRECISION RESULT >; END DFL ; DOUBLE PRECISION FIX FUNCTION "IDINT" ; DOUBLE TO INTEGER DEFINE IDF (X) < XALL ENTRY IDF.'X SIXBIT /IDF.'X/ IDF.'X: PUSH P,L ;SAVE THE SCRATCH REG HLRE L,X ;GET THE EXPONENT ASH L,-9 ;RIGHT 8 BITS JUMPGE X,IDF.XT ;JUMP IF POS. DMOVN X,X ;NEGATE TRC L,-1 ;COMPLEMENT THE EXPONENT IDF.XT: TLZ X,777000 ;CLEAR THE EXPONENT ASHC X,-201-^D26(L) ;CHANGE FRACTION TO INTEGER TLNE L,400000 ;SKIP IF POS. MOVN X,X ;NEGATE POP P,L ;RESTORE THE SCRATCH REG POPJ P, ;RETURN X=FIXED NUMBER >; END IDF ; DOUBLE PRECISION TO SINGLE FUNCTION DEFINE SNG (X)< XALL ENTRY SNG.'X SIXBIT /SNG.'X/ SNG.'X: JUMPL X,SNG3 ;NEGATIVE ARGUMENT? TLNE X+1,(1B1) ;POSITIVE. ROUND REQUIRED? TRON X,1 ;YES, TRY TO ROUND BY SETTING LSB POPJ P, ;WE WON, FINISHED MOVE X+1,X ;COPY HIGH PART OF ARG AND X,[777000,,1] ;MAKE UNNORMALIZED LSB, SAME EXPONENT FAD X,X+1 ;ROUND & RENORMALIZE POPJ P, ;HERE IF ARG IS NEGATIVE SNG3: DMOVN X,X ;MAKE POSITIVE TLNE X+1,(1B1) ;NEED ROUNDING? TRON X,1 ;YES, TRY TO DO IT BY SETTING LSB JRST SNG4 ;DONE MOVN X+1,X ;MAKE RE-NEGATED COPY OF HIGH PART ORCA X,[777,,-1] ;GET UNNORM NEG LSB WITH SAME EXPONENT FADR X,X+1 ;ROUND & NORMALIZE POPJ P, SNG4: MOVN X,X ;RE-NEGATE POPJ P, ;EXIT >; END SNG >; END IF1 PRGEND TITLE DFL.0 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 DFL 0 PRGEND TITLE DFL.2 SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 DFL 2 PRGEND TITLE DFL.3 SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 DFL 3 PRGEND TITLE DFL.4 SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 DFL 4 PRGEND TITLE DFL.6 SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 DFL 6 PRGEND TITLE DFL.10 SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 DFL 10 PRGEND TITLE DFL.12 SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 DFL 12 PRGEND TITLE DFL.14 SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 DFL 14 PRGEND TITLE IDF.0 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 IDF 0 PRGEND TITLE IDF.2 SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 IDF 2 PRGEND TITLE IDF.4 SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 IDF 4 PRGEND TITLE IDF.6 SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 IDF 6 PRGEND TITLE IDF.10 SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 IDF 10 PRGEND TITLE IDF.12 SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 IDF 12 PRGEND TITLE IDF.14 SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 IDF 14 PRGEND TITLE SNG.0 ;COPYRIGHT (C) 1979,1981 BY DIGITAL EQUIPMENT CORPORATION ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 SNG 0 PRGEND TITLE SNG.2 SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 SNG 2 PRGEND TITLE SNG.4 SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 SNG 4 PRGEND TITLE SNG.6 SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 SNG 6 PRGEND TITLE SNG.10 SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 SNG 10 PRGEND TITLE SNG.12 SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 SNG 12 PRGEND TITLE SNG.14 SEARCH FORPRM,FORDAR NOSYM TWOSEG 400000 SNG 14 PRGEND TITLE DDIM Double precision positive difference ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED ; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE. ; ;COPYRIGHT (C) 1981 BY DIGITAL EQUIPMENT CORPORATION, MAYNARD, MASS. SEARCH FORPRM TWOSEG 400000 NOSYM ENTRY DDIM EXTERN DDIM. DDIM=DDIM. PRGEND TITLE DDIM. Double Precision Positive Difference. SUBTTL C. McCutcheon 6/26/81 ;COPYRIGHT (C) 1981 BY DIGITAL EQUIPMENT CORPORATION, MAYNARD, MASS. ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. ; Fortran Intrinsic Function for Positive Difference. ; Passed and returned arguments are double precision. ; Passed: A1,A2 ; Returns: ; A1-A2 if A1 .GT. A2 ; 0 if A1 .LE. A2 SEARCH FORPRM TWOSEG 400000 SALL HELLO (DDIM,.) PUSH P,T2 PUSH P,T3 DMOVE T0,@0(L) ;A1 DMOVE T2,@1(L) ;A2 CAMN T0,T2 ;If A1 .LT. A2, CAML T1,T3 CAMGE T0,T2 JRST RET0 ; return 0. DFSB T0,T2 ;A1-A2 RET: POP P,T3 POP P,T2 GOODBYE ;Return RET0: SETZB T0,T1 ;Return zero JRST RET PRGEND TITLE DINT Double precision truncation ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED ; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE. ; ;COPYRIGHT (C) 1981 BY DIGITAL EQUIPMENT CORPORATION, MAYNARD, MASS. SEARCH FORPRM TWOSEG 400000 NOSYM ENTRY DINT EXTERN DINT. DINT=DINT. PRGEND TITLE DINT. Double Precision Truncation. SUBTTL C. McCutcheon - 6/25/81 ;COPYRIGHT (C) 1981 BY DIGITAL EQUIPMENT CORPORATION, MAYNARD, MASS. ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. ; Fortran Instrinsic Function to return a Double Precision ; truncation of the Double Precision number passed. ; If Magnitude(A) .LT. 1, then 0 is returned, ; otherwise return the largest integer that does not exceed the ; magnitude of A, and whose sign is the same as that of A. SEARCH FORPRM TWOSEG 400000 SALL HELLO (DINT,.) ;Enter routine DMOVE T0,@0(L) ;Get argument to truncate. JUMPN T0,NZERO ;If zero GOODBYE ;Return. NZERO: PUSH P,T2 ;Save ac's SKIPGE @0(L) ;If original number .LT. zero, DMOVN T0,T0 ; negate it. HLRZ T2,T0 ;Get exponent LSH T2,-^D9 ;Put rightmost CAIG T2,200 ;If number .LT. 1.0D0, JRST ZERO ; return zero. CAIL T2,276 ;If exponent .GE. 276 then there is JRST DONE ; no fraction, return the number passed. ; Now shift out the fractional part of the number. ASHC T0,-276(T2) ;Shift into integer position. MOVN T2,T2 ;Negate exponent. ASHC T0,276(T2) ;Shift back to where found. SKIPGE @0(L) ;If original number .LT. zero, DMOVN T0,T0 ; negate it again. BYE: POP P,T2 ;Pop saved ac's GOODBYE ;Return ; No fractional part, return number passed to function DONE: DMOVE T0,@0(L) JRST BYE ZERO: SETZB T0,T1 ;Return a 0. JRST BYE PRGEND TITLE DNINT Double precision nearest whole number ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED ; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE. ; ;COPYRIGHT (C) 1981 BY DIGITAL EQUIPMENT CORPORATION, MAYNARD, MASS. SEARCH FORPRM TWOSEG 400000 NOSYM ENTRY DNINT EXTERN DNINT. DNINT=DNINT. PRGEND TITLE DNINT. Double Precision nearest whole number. SUBTTL C. McCutcheon - 6/25/81 ;COPYRIGHT (C) 1981 BY DIGITAL EQUIPMENT CORPORATION, MAYNARD, MASS. ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. ; Fortran Instrinsic Function to return the double precision ; nearest whole number to the double precision number passed. ; Number returned is: Integer(A+.5) if A .GE. 0 ; Integer(A-.5) if A .LT. 0. SEARCH FORPRM TWOSEG 400000 HELLO (DNINT,.) ;Enter routine DMOVE T0,@0(L) ;Get argument to truncate. JUMPN T0,NZERO ;Is number passed = 0? GOODBYE ;Yes. NZERO: PUSH P,T2 ;Save an ac CAIGE T0,0 ;If original number .LT. zero, DMOVN T0,T0 ; negate it. TLNN T0,200000 ;Is number .LT. 1/2? JRST ZERO ;Yes, go return zero. HLRZ T2,T0 ;Get exponent LSH T2,-^D9 ;Put rightmost CAIL T2,276 ;If no fraction possible, JRST DONE ; return the number passed. ; Now shift out the fractional part of the number. DFAD T0,[200400,,0 0,,0] ;Add .5D0 before truncation. HLRZ T2,T0 ;Get new exponent. LSH T2,-^D9 ;Put rightmost. ASHC T0,-276(T2) ;Shift into integer position. MOVN T2,T2 ;Negate exponent. ASHC T0,276(T2) ;Shift back to where found. SKIPGE @0(L) ;If original number .LT. zero, DMOVN T0,T0 ; negate it again. BYE: POP P,T2 ;Pop saved ac. GOODBYE ;Return. ; No fractional part found, return number passed to function. DONE: DMOVE T0,@0(L) JRST BYE ZERO: SETZB T0,T1 ;Return 0. JRST BYE PRGEND TITLE DPROD Double precision product for single precision ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED ; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE. ; ;COPYRIGHT (C) 1981 BY DIGITAL EQUIPMENT CORPORATION, MAYNARD, MASS. SEARCH FORPRM TWOSEG 400000 NOSYM ENTRY DPROD EXTERN DPROD. DPROD=DPROD. PRGEND TITLE DPROD. Double Precision product for single precision. SUBTTL C. McCutcheon - 6/25/81 ;COPYRIGHT (C) 1981 BY DIGITAL EQUIPMENT CORPORATION, MAYNARD, MASS. ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. ; Fortran Instrinsic Function to return a double precision ; product for two real arguments passed it. SEARCH FORPRM TWOSEG 400000 SALL HELLO (DPROD,.) ;Enter routine MOVE T1,@1(L) ;Get 2nd argument MOVEM T1,MULT ;Save 2nd argument MOVE T0,@0(L) ;Get 1st argument SETZB T1,MULT+1 ;Set both second words to 0 DFMP T0,MULT ;Multiply and leave result in AC0. GOODBYE ;Return RELOC MULT: BLOCK 2 ;Multiplier RELOC PRGEND TITLE IDNINT Integer nearest whole number for double precision ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED ; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE. ; ;COPYRIGHT (C) 1981 BY DIGITAL EQUIPMENT CORPORATION, MAYNARD, MASS. SEARCH FORPRM TWOSEG 400000 NOSYM ENTRY IDNINT EXTERN IDNIN. IDNINT=IDNIN. PRGEND TITLE IDNIN. Integer nearest whole number for double prec. SUBTTL C. McCutcheon - 6/25/81 ;COPYRIGHT (C) 1981 BY DIGITAL EQUIPMENT CORPORATION, MAYNARD, MASS. ;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY BE USED AND COPIED ;ONLY IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE AND WITH THE ;INCLUSION OF THE ABOVE COPYRIGHT NOTICE. THIS SOFTWARE OR ANY OTHER ;COPIES THEREOF MAY NOT BE PROVIDED OR OTHERWISE MADE AVAILABLE TO ANY ;OTHER PERSON. NO TITLE TO AND OWNERSHIP OF THE SOFTWARE IS HEREBY ;TRANSFERRED. ;THE INFORMATION IN THIS SOFTWARE IS SUBJECT TO CHANGE WITHOUT NOTICE ;AND SHOULD NOT BE CONSTRUED AS A COMMITMENT BY DIGITAL EQUIPMENT ;CORPORATION. ;DIGITAL ASSUMES NO RESPONSIBILITY FOR THE USE OR RELIABILITY OF ITS ;SOFTWARE ON EQUIPMENT WHICH IS NOT SUPPLIED BY DIGITAL. ; Fortran Intrinsic Function to return the nearest ; integer to the double precision number passed. ; Number returned is: Integer(A+.5) if A .GE. 0 ; Integer(A-.5) if A .LT. 0. SEARCH FORPRM TWOSEG 400000 SALL HELLO (IDNIN.,IDNINT) ;Enter routine DMOVE T0,@0(L) ;Get argument to round. JUMPN T0,NZERO ;If number passed is 0, GOODBYE ; return NZERO: CAIGE T0,0 ;If original number .LT. zero, DMOVN T0,T0 ; negate it. TLNN T0,200000 ;Is |number| .LT. 1/2 ? JRST ZERO ;Yes, go return zero PUSH P,T2 ;Save an ac DFAD T0,[200400,,0 ;Add 0.5D0 before truncation 0,,0] ; to round the number ; Now shift out the fractional part of the number. HLRZ T2,T0 ;Get exponent. LSH T2,-^D9 ;Put rightmost. CAILE T2,243 ;If not representable, (overflow imminent) JRST DONE ; return largest integer TLZ T0,777000 ;Eliminate exponent. ASHC T0,-233(T2) ;Shift into integer position, leaving ; integer result in T0. SKIPGE @0(L) ;If original number .LT. zero, MOVN T0,T0 ; negate it a