Trailing-Edge
-
PDP-10 Archives
-
AP-D480B-SB_1978
-
fordbl.mac
There are 3 other files named fordbl.mac in the archive. Click here to see a list.
TITLE DABS %4.(235) DOUBLE PRECISION ABSOLUTE VALUE
SUBTTL H. P. WEISS/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1973,1977 BY DIGITAL EQUIPMENT CORPORATION
SEARCH FORPRM
ENTRY DABS
EXTERN DABS.
DABS=DABS.
PRGEND
TITLE DABS. %4.(235) DOUBLE PRECISION ABSOLUTE VALUE
SUBTTL D. TODD /KK/DMN/DRT/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1972,1977 BY DIGITAL EQUIPMENT CORPORATION
;DOUBLE PRECISION ABSOLUTE VALUE FUNCTION
;THIS ROUTINE RETURNS THE ABSOLUTE VALUE OF A DOUBLE PRECISION
;ARGUMENT
;THE CALLING SEQUENCE FOR THE ROUTINE IS
; JSA Q, DABS
; EXP ARG
;WHERE ARG IS THE ADDRESS OF THE HIGH ORDER PART OF A DOUBLE
;PRECISION ARGUMENT, THE LOW ORDER PART BEING IN ARG+1. THE
;DOUBLE PRECISION ANSWER IS RETURNED IN ACCUMULATORS A AND B.
SEARCH FORPRM
A= 0
B= 1
Q= 16
P= 17
HELLO (DABS,.) ;[235] ENTRY TO DABS ROUTINE
DMOVE A,@(Q) ;GET ARGUMENT
SKIPGE A ;IS ARGUMENT POSITIVE?
DFN A, B ;NO, NEGATE IT
GOODBY (1) ;EXIT
PRGEND
TITLE DMOD %4.(235) DOUBLE PRECISION MOD FUNCTION
SUBTTL H. P. WEISS/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1973,1977 BY DIGITAL EQUIPMENT CORPORATION
SEARCH FORPRM
ENTRY DMOD
EXTERN DMOD.
DMOD=DMOD.
PRGEND
TITLE DMOD. %4.(235) DOUBLE PRECISION MOD FUNCTION
SUBTTL D. TODD /DRT/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1972,1977 BY DIGITAL EQUIPMENT CORPORATION
;THIS ROUTINE CALCULATES
; MOD(A,B) = A-[A/B]*B
;FOR DOUBLE PRECISION ARGUMENTS A AND B, WHERE [A/B] IS THE
;GREATEST INTEGER IN THE MAGNITUDE OF A/B.
;THE CALLING SEQUENCE FOR THE ROUTINE IS AS FOLLOWS:
; JSA Q, DMOD
; EXP ARG1
; EXP ARG2
;ARG1 AND ARG2 ARE THE ADDRESSES OF THE HIGH ORDER PORTIONS OF
;THE DOUBLE PRECISION ARGUMENTS. THE DOUBLE PRECISION ANSWER IS
;RETURNED IN ACCUMULATORS A AND B.
SEARCH FORPRM
P=17 ;PUSH DOWN POINTER
Q=16 ;JSA-JRA ACCUMULATOR
HELLO (DMOD,.) ;[235] ENTRY TO DMOD ROUTINE.
DMOVE 0,@1(Q) ;GET ARG B
SKIPGE 0 ;GET
DFN 0,1 ;/ARG B/.
DMOVEM 0,ARGB ;SAVE ARG B IN ARGB.
DMOVE 0,@(Q) ;GET ARG A
SKIPGE 0 ;GET
DFN 0,1 ;/ARG A/.
DMOVEM 0,ARGA ;SAVE ARG A IN ARGA.
MOVEM 2,SAV2 ;SAVE AC 2.
IFE CPU-KA10,<FDVL 0,ARGB ;DIVIDE A BY B
JFCL OVUNFL ;IF
MOVN 2,0 ;OVER
FMPR 2,ARGB+1 ;OR
JFCL ;UNDER
UFA 1,2 ;FLOW
FDVR 2,ARGB ;OCCURS,
JFCL ;GO
FADL 0,2 ;TO
JFCL OVUNFL > ;OVUNFL.
IFE CPU-KI10,<DFDV 0,ARGB
JFCL OVUNFL >
DMOVEM 0,SAVNUM ;SAVE NUMBER
LDB 2,[POINT 8,0,8] ;GET EXPONENT
TRC 2,777777 ;GET 1'S COMPLEMENT OF EXPONENT
MOVSI 0,777000 ;INITIALIZE MASK
MOVEI 1,0 ;IN 0 AND 1.
CAIL 2,777577 ;IS THE NUMBER ALL FRACTION BITS?
TDZA 0,0 ;YES, A FRACTION WILL TRUNCATE TO 0
ASHC 0,201(2) ;CREATE MASK(SHIFT TO RIGHT ONLY)
;REMEMBER- EXPONENT IS 1'S COMP NEGATIVE
IFE CPU-KA10,<ASH 1,-8 > ;PROTECT LOW ORDER EXPONENT
AND 0,SAVNUM ;MASK HI WORD
AND 1,SAVNUM+1 ;AND LOW WORD
IFE CPU-KA10,<FADL 0,1 > ;AND STRAIGHTEN 0'S OUT.
FLMUL 0,ARGB ;CALCULATE [A/B]*B.
DFN 0,1 ;GET -[A/B]*B
FLADD 0,ARGA ;CALC. A-[A/B]*B.
MOVE 2,SAV2 ;RESTORE AC 2.
SKIPGE @(Q) ;GIVE THE ANSWER
DFN 0,1 ;THE CORRECT SIGN.
GOODBY (2) ;EXIT.
OVUNFL: MOVE 2,SAV2 ;RESTORE AC 2.
JUMPE 0,ANSISA ;IF UNDERFLOW, GO TO ANSISA.
ERROR (APR,5,1,.+1) ;OVERFLOW. RETURN AN
SETZB 0,1 ;AN ANSWER OF ZERO AND
GOODBY (2) ;EXIT.
ANSISA: DMOVE 0,@(Q) ;UNDERFLOW. ANS = A
GOODBY (2)
ARGA: BLOCK 2
ARGB: BLOCK 2
SAVNUM: BLOCK 2
SAV2: BLOCK 1
PRGEND
TITLE DMAX1 %4.(235) DOUBLE PRECISION MAXIMUM OF A SERIES OF ARGUMENTS
SUBTTL H. P. WEISS/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1973,1977 BY DIGITAL EQUIPMENT CORPORATION
SEARCH FORPRM
ENTRY DMAX1
EXTERN DMAX1.
DMAX1=DMAX1.
PRGEND
TITLE DMAX1. %4.(235) DOUBLE PRECISION MAXIMUM OF A SERIES OF ARGUMENTS
SUBTTL D. TODD /DRT/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1972,1977 BY DIGITAL EQUIPMENT CORPORATION
;THIS ROUTINE IS CALLED IN THE FOLLOWING MANNER
; JSA Q, DMAX1
; EXP ARG1
; EXP ARG2
; .
; .
; .
;WHERE ARG1,ARG2,...ARE THE ADDRESSES OF THE HIGH ORDER
;PORTIONS OF DOUBLE PRECISION ARGUMENTS. THE MAXIMUM OF THE
;ENTIRE SET IS RETURNED AS A DOUBLE PRECISION NUMBER IN
;ACCUMULATOR A AND B.
SEARCH FORPRM
A= 0
B= 1
C= 14
Q= 16
HELLO (DMAX1,.) ;[235] ENTRY TO DMAX1 ROUTINE
PUSH P,C ;SAVE AC C
IFN F40LIB,<
PUSH P,L ;SAVE THE LINK FOR F40
TLZN L,-1 ;F40 CALL
>
HLL L,-1(L) ;GET THE F10 ARG COUNT
DMOVE A,@(Q) ;GET FIRST ARGUMENT
JRST DMAX.2 ;ADDRESS OF NEXT, START CHECKING
DMAX.1: MOVEI C,@0(Q) ;GET ADDRESS OF NEXT ARGUMENT
CAMLE A,(C) ;IS HIGH ORDER > THAN THIS ONE?
JRST DMAX.2 ;YES, GET ADDRESS OF NEXT IN LIST
CAME A,(C) ;ARE HIGH ORDER WORDS EQUAL?
JRST DMAX.3 ;NO, REPLACEMENT NECESSARY
IFN CPU-KI10,<CAML B, 1(C) ;YES, CHECK LOW ORDER WORDS
JRST DMAX.2 > ;OK, GET ADDRESS OF NEXT ONE
IFE CPU-KI10,<CAMGE B,1(C) > ;SKIP IF PRESENT ARG IS LARGER
DMAX.3: DMOVE A,(C) ;PICK UP NEW ARG
DMAX.2: AOBJN L,DMAX.1 ;END OF ARG LIST
IFN F40LIB,<
TLNN L,-1 ;F40 CALL
JRST DMAX.4 ;NO, F10 CALL EXIT
MOVE C,(L) ;GET THE NEXT ARGUMENT
TLC C,(<JUMP>) ;COMPLEMENT OP-CODE "JUMP" BITS
TLNN C,777000 ;IS THE ARG "JUMP" ?
JRST DMAX.1 ;YES,INCREMENT AND CHECK FURTHER
DMAX.4:
POP P,L ;RESTORE LINK
>
POP P,C ;NO, RESTORE C
GOODBY (2)
PRGEND
TITLE DMIN1 %4.(235) DOUBLE PRECISION MINIMUM OF A SERIES OF ARGUMENTS
SUBTTL H. P. WEISS/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1973,1977 BY DIGITAL EQUIPMENT CORPORATION
SEARCH FORPRM
ENTRY DMIN1
EXTERN DMIN1.
DMIN1=DMIN1.
PRGEND
TITLE DMIN1. %4.(235) DOUBLE PRECISION MINIMUM OF A SERIES OF ARGUMENTS
SUBTTL D. TODD /DRT/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1972,1977 BY DIGITAL EQUIPMENT CORPORATION
;THIS ROUTINE IS CALLED IN THE FOLLOWING MANNER:
; JSA Q, DMIN1
; EXP ARG1
; EXP ARG2
; .
; .
; .
;WHERE ARG1,ARG2,... ARE THE ADDRESSES OF THE HIGH ORDER
;PORTIONS OF DOUBLE PRECISION ARGUMENTS. THE MINIMUM OF
;THE ENTIRE SET IS RETURNED AS A DOUBLE PRECISION NUMBER
;IN ACCUMULATOR A AND B
SEARCH FORPRM
LIBSEG ;GET THE SEGMENT CONTROL
A= 0
B= 1
C= 14
Q= 16
HELLO (DMIN1,.) ;[235] ENTRY TO DMIN1 ROUTINE
PUSH P,C ;SAVE AC C
IFN F40LIB,<
PUSH P,L ;SAVE THE LINK
TLZN L,-1 ;F40 CALL
>
HLL L,-1(L) ;NO F10 GET THE ARG COUNT
DMOVE A,@(Q) ;GET FIRST ARGUMENT
JRST DMIN.2 ;ADDRESS OF NEXT,START CHECKING
DMIN.1: MOVEI C,@0(Q) ;GET ADDRESS OF NEXT ARG
CAMGE A, (C) ;IS HIGH ORDER LESS THAN NEXT ARG?
JRST DMIN.2 ;YES, GET ADDRESS OF NEXT
CAME A, (C) ;NO, ARE HIGH ORDER WORDS EQUAL?
JRST DMIN.3 ;NO, REPLACEMENT IS NECESSARY
IFN CPU-KI10,<CAMG B, 1(C) ;YES, CHECK LOW ORDER WORDS
JRST DMIN.2 > ;PRESENT ARGUMENT IS SMALLER
IFE CPU-KI10,<CAMLE B,1(C) > ;SKIP IF PRESENT ARG IS SMALLER
DMIN.3: DMOVE A,(C) ;PICK UP NEW ARG
DMIN.2: AOBJN L,DMIN.1 ;CONTINUE THRU THE LIST
IFN F40LIB,<
TLNN L,-1 ;F40 CALL
JRST DMIN.4 ;NO, F10 CALL EXIT
MOVE C,(L) ;GET THE NEXT ARGUMENT
TLC C,(<JUMP>) ;COMPLEMENT OP CODE "JUMP" BITS
TLNN C,777000 ;IS THE ARG "JUMP"?
JRST DMIN.1 ;YES , INCREMENT AND CHECK FURTHER
DMIN.4:
POP P,L ;RESTORE THE LINK
>
POP P,C ;NO, RESTORE AC C AND EXIT
GOODBY (2)
PRGEND
TITLE DSIGN %4.(235) DOUBLE PRECISION SIGN ROUTINE
SUBTTL H. P. WEISS/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1973,1977 BY DIGITAL EQUIPMENT CORPORATION
SEARCH FORPRM
ENTRY DSIGN
EXTERN DSIGN.
DSIGN=DSIGN.
PRGEND
TITLE DSIGN. %4.(235) DOUBLE PRECISION SIGN ROUTINE
SUBTTL D. TODD /DRT/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1972,1977 BY DIGITAL EQUIPMENT CORPORATION
;FROM 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
;THE CALLING SEQUENCE FOR THIS ROUTINE IS AS FOLLOWS
; JSA Q,DSIGN
; EXP ARG1
; EXP ARG2
;ARG1 AND ARG2 ARE THE ADDRESSES OF THE HIGH ORDER WORDS OF
;THE DOUBLE PRECISION ARGUMENTS, THE DOUBLE PRECISION
;ANSWER IS RETURNED IN ACCUMULATORS A AND B.
SEARCH FORPRM
A=0
B=1
Q=16
HELLO (DSIGN,.) ;[235] ENTRY TO DSIGN ROUTINE
DMOVE A,@(Q) ;GET FIRST ARGUMENT
SKIPGE @1(Q) ;THEN
JUMPL A,OUT ;CHOOSE
SKIPL @1(Q) ;THE
JUMPGE A,OUT ;CORRECT
DFN A,B ;SIGN.
OUT: GOODBY (2) ;EXIT
PRGEND
TITLE DEXP.3 %5A(643) PDP-10/10I DOUBLE PRECISION EXP.3 FUNCTION
SUBTTL D. TODD /DRT/SWG 24-FEB-1977
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1972,1977 BY DIGITAL EQUIPMENT CORPORATION
;THIS PROGRAM CALCULATES A**B, WHERE A AND B ARE DOUBLE
;PRECISION FLOATING POINT NUMBERS. THE ALGORITHM USED IS
;A**B = DEXP(B*DLOG(/A/)). THE ABSOLUTE VALUE OF A IS USED
;IN THIS CALCULATION BECAUSE A NEGATIVE NUMBER TO A
;NON-INTEGER POWER PRODUCES A COMPLEX ANSWER, AND B IS
;PRESUMED TO BE NON-INTEGER.
;[624]CHECKS TO SEE IF B IS AN INTEGER AND TREATS IT
;[624] ACCORDINGLY BY CALLING DEXP.2. GIVES LIBRARY ERROR FOR
;[624]NEG NUM TO NON-INTEGER POWER.
;THE CALLING SEQUENCE FOR THE ROUTINE IS AS FOLLOWS:
; MOVEI Q, ARG2
; PUSHJ P, DEXP3
;ARG2 IS THE ADDRESS OF THE HIGH ORDER PART OF THE DOUBLE
;PRECISION POWER, AND THE BASE, ARG1, IS IN ACCUMULATORS
;A AND B WHEN THE ROUTINE IS CALLED. THE DOUBLE PRECISION
;ANSWER IS LEFT IN ACCUMULATORS A AND B.
;** ; [624] INSERT BEFORE SEARCH - SWG 29 NOV 1976
EXTERN DEX2.. ;[624]STANDARD DP TO INTEGER POWER
SEARCH FORPRM
IFN F40LIB,<ENTRY DEXP.3>
IFN F10LIB,<ENTRY DEXP3.>
E= 4 ;[624]
D= 3 ;[624]
C= 2
B= 1
A= 0
Q= 16
P= 17
;**; [643] DEXP3.-2 INSERT DEF OF ARGAX SWG 24-FEB-77
ARGAX: BLOCK 2 ;[643]
IFN F10LIB,<
SIXBIT /DEXP3./
DEXP3.: DMOVE T0,@(L) ;GET THE BASE
MOVEI L,@1(L)
IFN F40LIB,<
JRST DEXP.3 ;COMMON ROUTINE
>>
IFN F40LIB,<
SIXBIT /DEXP.3/
DEXP.3:>
;**;[665] INSERT @ DEXP.3+1L FROM BELOW SWG 28-JUL-77
PUSH P,C ;[665] SAVE TEMPORARY REGISTERS
PUSH P,D ;[665] AT ENTRY BEFORE DOING
PUSH P,E ;[665] ANYTHING FOR CLEAN EXIT
SKIPN (Q) ;IS EXPONENT ZERO?
JRST DEXPZ ;YES, RETURN ANSWER OF 1.0
JUMPE A,[SKIPL (Q) ;IS BASE 0 WITH POS EXPONENT?
;**; [665] CHANGE @ DEXP.3+4L SWG 28-JUL-77
JRST DRET ;[665]YES, RETURN 0
JRST OV4] ;NO, BASE 0 AND NEG EXP- OVERFLOW
;**; [665] MOVE @ DEXP.3+6L 3 LINES ABOVE SWG 28-JUL 77
;**; [624] DEXP.3+6L DELETE 2 LINES AND INSERT SWG 29-NOV-76
DMOVEM A,ARGAX ;[624]SAVE BASE
DMOVE A,(Q) ;[624]PULL IN EXPONENT
MOVM D,A ;[624]MAGNITUDE OF WORD 1 OF EXP
MOVEI C,0 ;[624]CLEAR AC C
LSHC C,11 ;[624]SHIFT 9 BITS
CAILE C,233 ;[624]INTEGER MUST FIT IN 27 BITS
JRST DEXPD ;[624]NOT AN INTEGER
SUBI C,200 ;[624]200 COMPLEMENT ON EXP
MOVM E,B ;[624]WORD 2 OF EXP
IFE CPU-KA10,< ;[624]
HRLZI B,777000 ;[624]WANT TO LOOK AT FRACTION
ANDCAM B,E ;[624]CLEAR FIRST 9 BITS
> ;[624]
JUMPN E,DEXPD ;[624]IS WORD 2 ZERO?
LSH D,-1 ;[624]COULD STILL BE INTEGER
HRRZ E,C ;[624]SET E AS INDEX REG
MOVEI C,0 ;[624]CLEAR C
ASHC C,(E) ;[624]SHIFT LEFT BY CONTENT S OF E
JFCL DEXPD ;[624]OVERFLOW
JUMPN D,DEXPD ;[624]IS IT INTEGER??
MOVE B,(Q) ;[624]YES-NEED SIGN OF EX
SKIPGE ,B ;[624]NEGATIVE??
MOVNS C ;[624]YES-NEGATE INTEGER
PUSH P,(Q) ;[624]SAVE WHAT'S POINTED TO BY Q
MOVEM C,(Q) ;[624]SET ARG FOR DEXP.2
DMOVE A,ARGAX ;[624]RESTORE BASE
PUSHJ P,DEX2.. ;[624]
POP P,(Q) ;[624]RESTORE ARG
JRST DRET ;[624]RESTORE REGS AND RETURN
DEXPD: DMOVE A,ARGAX ;[624]CHECK SIGN
JUMPGE A,DEXPD2 ;[624]NEGATIVE??
DFN A,B ;[624]YES-GET ABS VALUE
ERROR (LIB,2,2,[ASCIZ/Attempt to raise negative double precision number to non-integer power/]) ;[624]
DEXPD2: DMOVEM A,ARGAX ;[624-INSERT LABEL]SAVE BASE
FUNCT DLOG.,<ARGAX> ;CALC.
;LOG(/A/).
DMOVEM A,ARGAX ;STORE IT IN ARGAX.
DMOVE A,(Q) ;ARG B TO AC'S 0 AND 1.
MOVEM C,CSAVE ;SAVE AC 2.
IFE CPU-KA10,<MOVEM A,C ;CALCULATE
FMPR C,ARGAX+1 ;B*
JFCL ;LOG(/A/)
FMPR B,ARGAX ;AND
JFCL ;GO
UFA B,C ;TO
JFCL ;OVUNFL
FMPL A,ARGAX ;IF
JFCL OVUNFL ;OVER-
UFA B,C ;OR
FADL A,C ;UNDERFLOW
JFCL OVUNFL > ;OCCUR.
IFE CPU-KI10,<DFMP A,ARGAX
JFCL OVUNFL >
DMOVEM A,ARGAX ;STORE B*LOG(/A/) IN ARGAX.
FUNCT DEXP.,<ARGAX> ;CALC. EXP[B*LOG(/A/)] AND
;LEAVE IT IN AC'S 0 AND 1.
MOVE C,CSAVE ;RESTORE AC 2.
DRET: POP P,E ;[624]RESTORE REGS
POP P,D ;[624]
POP P,C ;[624]
POPJ P, ;EXIT.
DEXPZ: MOVSI A, (1.0) ;ANS = 1.0
MOVEI B, 0 ;AND
;**; [665] CHANGE @ DEXPZ+2 SWG 28-JUL-77
JRST DRET ;[665]EXIT.
OVUNFL: MOVE C,CSAVE ;RESTORE AC 2.
JUMPE A,DEXPZ ;IF EXP= 0, ANS = 1.0.
JUMPL A,OV6 ;GO TO UNDERFLOW.
OV4: ERROR (APR,5,1,.+1) ;OVERFLOW. OUTPUT
HRLOI A,377777 ;ANS =
IFE CPU-KA10,<HRLOI B,344777 > ;INFINITY.
IFN CPU-KA10,<HRLOI B,377777 >
;**; [665] CHANGE @ OV6-1 SWG 28-JUL-77
JRST DRET ;[665]EXIT.
OV6: ERROR (APR,7,1,.+1) ;UNDERFLOW. OUTPUT
SETZB A,B ;ANS = 0, AND
;**; [665] CHANGE @ OV6+2 SWG 28-JUL-77
JRST DRET ;[665]EXIT.
;**;[643] CSAVE-1 MOVE DEF OF ARGAX ABOVE SWG 24-FEB-77
CSAVE: BLOCK 1
PRGEND
TITLE DEXP.2 %5A(624) PDP-10/10I DOUBLE PRECISION EXP2 FUNCTION
SUBTTL D. TODD /DRT/SWG 29-NOV-1976
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1972,1977 BY DIGITAL EQUIPMENT CORPORATION
;THIS ROUTINE CALCULATES A DOUBLE PRECISION NUMBER RAISED
;TO A FIXED POINT POWER. THE CALCULATION IS A**N, WHERE N
;IS AN INTEGER OF THE FORM
; N= Q(0) + Q(1)*2 + Q(2)*4 + ... WHERE Q(I) = 0 OR 1
;THE ONLY RESTRICTION ON THE BASE OR EXPONENT IS THAT
;AN EXPONENT OF 400000000000 IS NOT HANDLED CORRECTLY.
;THE ROUTINE IS CALLED BY
; MOVEI Q, POWER
; PUSHJ P, DEXP2
;WHERE POWER IS THE ADDRESS OF THE FIXED POINT POWER, AND
;THE DOUBLE PRECISION BASE IS IN ACCUMULATORS A AND B. THE
;RESULT IS RETURNED IN ACCUMULATORS A AND B.
;**;[624] INSERT BEFORE SEARCH SWG 29-NOV-1976
ENTRY DEX2.. ;[624]ENTRY FROM DEXP.3
SEARCH FORPRM
IFN F40LIB,<ENTRY DEXP.2>
IFN F10LIB,<ENTRY DEXP2.>
A= 0
B= 1
C= 2
D= 3
E= 4
G= 6
Q= 16
P= 17
X= G ;HIGHEST AC TO SAVE
IFN F10LIB,<
SIXBIT /DEXP2./
DEXP2.: DMOVE T0,@(L) ;GET THE BASE
MOVEI L,@1(L) ;POINT TO THE POWER
IFN F40LIB,<
JRST DEXP.2 ;GO TO COMMON ROUTINE
>>
IFN F40LIB,<
SIXBIT /DEXP.2/
DEXP.2:>
;**; [624] INSERT LABEL AT DEXP.2+1 SWG 29-NOV-1976
DEX2..: SKIPN (Q) ;[624]IS EXPONENT 0?
JRST [MOVSI A,(1.0) ;YES, A**0 GIVES 1
MOVEI B,0
POPJ P,]
JUMPE A,[SKIPL (Q) ;IS BASE 0 WITH POSITIVE EXPONENT?
POPJ P, ;YES, RETURN 0
ERROR (APR,5,1,.+1) ;BASE IS 0 WITH NEG. EXP- OVERFLOW
HRLOI A,377777 ;RETURN LARGEST POSITIVE NUMBER
IFE CPU-KA10,<HRLOI B,344777 >
IFN CPU-KA10,<HRLOI B,377777 >
POPJ P,]
MOVEM X,XSAVE ;SAVE AC TO DO BLT
MOVE X,XBLT ;SAVE OTHER AC'S
BLT X,XSAVE-1 ;...
SKIPL G,(Q) ;GET EXPONENT. IS IT NEGATIVE?
JRST [DMOVE D,A ;NO, PUT ARG IN D,D+1
JRST DEX2] ;START MAIN LOOP
MOVMS G ;GET POSITIVE VALUE
MOVSI D,(1.0) ;GET DOUB. PRECISION 1.0
MOVEI E,0 ;...
;CALCULATE (1/X)**N, SINCE N .L. 0
IFE CPU-KA10,<FDVL 3,0
JFCL 1,OVER
MOVN 5,3
FMPR 5,1
JFCL
UFA 4,5
FDVR 5,0
JFCL
FADL 3,5 >
IFE CPU-KI10,<DFDV D,0
JFCL 1,OVER >
DEX2: MOVSI A,(1.0) ;GET DOUB. PREC. 1.0
MOVEI B,0 ;...
JRST DEX4 ;START CALCULATING POWERS OF X (OR 1/X)
DEX3: ;SQUARE X (OR 1/X) AGAIN
IFE CPU-KA10,<DMOVEM D,TEMP
MOVEM D,D+2
FMPR D+2,TEMP+1
JFCL
FMPR D+1,TEMP
JFCL
UFA D+1,D+2
JFCL
FMPL D,TEMP
JOV OVR
UFA D+1,D+2
FADL D,D+2
JOV OVR >
IFE CPU-KI10,<DFMP D,D
JOV OVR >
LSH G,-1 ;LOOK AT NEXT BIT IN N
DEX4: TRZN G,1 ;IS LO BIT IN N A 1?
JRST DEX5 ;NO, DON'T MULTIPLY INTO ANSWER
;MULTIPLY POWER OF X INTO ANSWER
IFE CPU-KA10,<MOVEM A,A+2
FMPR A+2,D+1
JFCL
FMPR A+1,D
JFCL
UFA A+1,A+2
JFCL
FMPL A,D
JOV DEX6
UFA A+1,A+2
FADL A,A+2
JOV DEX6 >
IFE CPU-KI10,<DFMP A,D
JOV DEX6 >
DEX5: JUMPN G,DEX3 ;IF G .N. 0, GET MORE POWERS OF X (OR 1/X)
DEX6: MOVS X,XBLT ;RESTORE AC'S
BLT X,X ;...
CPOPJ: POPJ P,
OVR: ;ARITHMETIC FAULT, MOVE FIX UP TO A,B
SKIPGE A ;SHOULD RESULT BE NEGATIVE?
DFN D,E ;YES
OVR2:
DMOVE A,D
JRST DEX6 ;AND EXIT
OVER: JUMPGE D,OVR2 ;IF THE ARG IS <0 AND THE EXPONENT
TRNN G,1 ;IS ODD, THEN
DFN D,E ;THE ANSWER
JRST OVR2 ;IS < 0.
XBLT: XWD C,ACSAVE ;BLT POINTER TO RESTORE AC'S
ACSAVE: BLOCK X-C
XSAVE: BLOCK 1 ;FOR AC X
TEMP: BLOCK 2
PRGEND
TITLE DEXP %4.(235) PDP-10/10I DOUBLE PRECISION EXPONENTIAL FUNCTION
SUBTTL H. P. WEISS/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1973,1977 BY DIGITAL EQUIPMENT CORPORATION
SEARCH FORPRM
ENTRY DEXP
EXTERN DEXP.
DEXP=DEXP.
PRGEND
TITLE DEXP. %5A(700) PDP-10/10I DOUBLE PRECISION EXPONENTIAL FUNCTION
SUBTTL D. TODD /DRT/HPW/MD 12-AUG-74
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1972,1977 BY DIGITAL EQUIPMENT CORPORATION
;ARGUMENT. AN ARGUMENT OF ZERO CAUSES AN IMMEDIATE EXIT WITH AN
;ANSWER OF 1.0 . AN ARGUMENT WHOSE MAGNITUDE EXCEEDS 88.028
;CAUSES THE ROUTINE TO EXIT WITH 0 IF THE ARGUMENT WAS NEGATIVE.
;AND 377777777777 IF THE ARGUMENT WAS POSITIVE. THIS IS
;BECAUSE DIRECT CALCULATION OF EXP(X) FOR ABSF(X)>88.028 WOULD
;CAUSE EXPONENT OVERFLOW OR UNDERFLOW.
;THE ROUTINE USES THE FOLLOWING ALGORITHM:
;EXP(X) = 2**(X*LOG2(E))
; = 2**(M+F) WHERE M IS AN INTEGER AND 0<F<1
; = 2**(M+N+R) WHERE 0<R<1/8 AND M+N+R=X*LOG2(E)
; = 2**(M+N) * EXP(R*LOG(2))
;2**M IS CALCULATED EASILY WITH THE FLOATING SCALE INSTRUCTION.
;2**N IS CALCULATED BY DETERMINING THE CORRECT INTERVAL OF N AND
;USING A TABLE OF POWERS OF TWO FROM 2**1/8 TO 2**7/8.
;FINALLY, EXP(R*LOG(2)) IS CALCULATED BY A CONTINUED FRACTION
;TAKEN FROM RALSTON AND WILF, "METHODS FOR DIGITAL COMPUTERS" :
;EXP(R*LOG(2)) = 1+A4/((B4/R) -C4 + D4*R + H4/(R + B4/R))
;THE FOLLOWING ERRORS HAVE BEEN OBSERVED WITH DEXP:
; 1. WITH -10.0<X<10.0, ERRORS RANGED FROM 0 TO 48
; UNITS IN THE 19TH SIGNIFICANT DIGIT. THE MEAN ERROR
; FOR 20 READINGS WAS 5.4 UNITS IN THE 19TH DIGIT.
; 2. WITH 10.0<X<40.0, ERRORS RANGED FROM 0 TO 41 UNITS
; IN THE 19TH SIGNIFICANT DIGIT. THE MEAN ERROR FOR
; 7 READINGS WAS 16 UNITS IN THE 19TH SIGNIFICANT DIGIT.
; 3. WITH 40.0<X<88.0, ERRORS RANGED FROM 5 TO 57 UNITS IN
; THE 19TH SIGNIFICANT DIGIT. THE MEAN ERROR FOR 12
; READINGS WAS 44.4 UNITS IN THE 19TH SIGNIFICANT DIGIT.
;THE ERRORS REFERRED TO ABOVE ARE ABSOLUTE ERRORS. IT SHOULD
;BE NOTED THAT ADDITIONAL ERRORS ARE INTRODUCED BY ERRORS IN
;THE DOUBLE PRECISION INPUT AND OUTPUT ROUTINES.
;700 23542 FIX FLOATING DIVIDE CHECK WHEN LEFT WORD
; OF MANTISSA IS <377 (OCTAL)
;****** END OF REVISION HISTORY
A= 0
B= 1
C= 2
D= 3
E= 4
F= 5
G= 6
Q= 16
P= 17
X= G ;HIGHEST AC TO SAVE
SEARCH FORPRM
;CONSTANTS AND TEMPORARY LOCATIONS AND STUFF
XBLT: XWD C,ACSAVE
DCON1: 88.028
DLOG2E: DOUBLE 201561250731,112701376057 ;LOG2(E) = 1.44269 50408 88963 40740
TABLE: 0 ;0
040000000000 ;1/8
100000000000 ;2/8
140000000000 ;3/8
200000000000 ;4/8
240000000000 ;5/8
300000000000 ;6/8
340000000000 ;7/8
ONE: ;DOUBLE PRECISION 1.0
POWERS: DOUBLE 1.0,0 ;2**0 = 1.0
DOUBLE 201427127017,037250572672 ;2**1/8 = 1.09050 77326 65257 65919
DOUBLE 201460337602,214333425134 ;2**2/8 = 1.18920 71150 02721 06671
DOUBLE 201513773265,115425047073 ;2**3/8 = 1.29683 95546 51009 66590
DOUBLE 201552023631,237635714441 ;2**4/8 = 1.41421 35623 73095 04878
DOUBLE 201612634520,212520333270 ;2**5/8 = 1.54221 08254 07940 824
DOUBLE 201656423746,126551655275 ;2**6/8 = 1.68179 28305 07429 086
DOUBLE 201725403067,076722207113 ;2**7/8 = 1.83400 80864 09342 463
A4: DOUBLE 206744575555,062215755376 ;A4 = 60.59319 17173 36463 11080
B4: DOUBLE 207535527021,213670572221 ;B4 = 87.41749 72022 35527 474
MC4: DOUBLE 572033202222,715562022402 ;MC4=-C4 = -30.29659 58586 68231 555
D4: DOUBLE 201414631463,063146314632 ;D4 = 1.05
H4: DOUBLE 210654261010,261565402456 ;H4 = 214.17286 81454 77042 3113
ACSAVE: BLOCK X-C+1
B4F: BLOCK 2 ;TEMP FOR B4*F
FB4F: BLOCK 2 ;TEMP FOR F+B4*F
HELLO (DEXP,.) ;[235] ENTRY TO DEXP ROUTINE
MOVE 0,XBLT ;SAVE ACCUMULATORS
BLT 0,ACSAVE+X-C ;...
DMOVE A,@(Q) ;PICK UP ARGUMENT
JUMPE A,[MOVSI A,(1.0) ;RETURN 1.0 FOR ARG OF ZERO
JRST DEXEND] ;EXIT
MOVM C,A ;GET POS VALUE OF EXPONENT
CAML C,DCON1 ;TOO BIG TO COMPUTE?
JRST [HRLOI A,377777 ;YES, GET LARGEST POS NUM
MOVE C,ACSAVE ;RESTORE AC C.
SKIPGE @(Q) ;SUPPOSED TO BE SMALL ?
MOVEI A,1 ;YES, MAKE IT VERY SMALL
FADL A,A ;CAUSE OVER/UNDERFLOW, RETURN FIX UPS
JRST DEXEND] ;RETURN
FLMUL A,DLOG2E
HLRE E,A ;EXTRACT EXPONENT
ASH E,-9 ;...
TSC E,E ;TAKE 1'S COMPLEMENT IF NUM .L. 0
IFE CPU-KA10,<LSH B, 8 > ;REMOVE LOW ORDER EXP.
SKIPGE A ;CHANGE EXPONENT BITS TO SIGN BITS
TLOA A,377000 ;NUMBER NEGATIVE, SET BITS
TLZ A,377000 ;NUMBER POSITIVE, CLEAR BITS
ASHC A, 8 ;LEFT JUSTIFY ARG FRACTION BITS
;GET ANOTHER COPY OF FRACTION
DMOVE C,A
ASHC A, -243(E) ;SIMULATE THREE-WORD SHIFT WITH...
;TWO LONG SHIFTS. THIS LEAVES INTEGER
;IN A, FRAC IN B AND C.
LSH D,1 ;SQUEEZE OUT SIGN BIT
LSHC C,43-200(E) ;THEN DO 2ND LONG SHIFT, (THE LSHC HERE
;PREVENTS OVERFLOW GOING LEFT)
TLZ B, (1B0) ;CLEAR SIGN BIT. IF FRAC WAS <0,
;THIS GIVES 1-FRAC AND PROPER EXP.
HRRM A, SCALE ;SAVE EXPONENT FOR FUTURE SCALING
MOVEI G, 7 ;GET INDEX REGISTER POINTER TO TABLE
;**; [700] CHANGE & INSERT @ REDUCE SWG 30-AUG-77
REDUCE: CAME B, TABLE(G) ;[700]DOES ARGUMENT MATCH TABLE ENTRY?
JRST REDUC2 ;[700] NO - KEEP LOOKING
JFFO C,.+2 ;[700] HOW MANY LEADING 0'S
JRST REDUC1 ;[700] NONE - C CONTAINS 0
CAIL D,^D28 ;[700] ARE ANY LEFT OF BIT 28 ON?
REDUC1: JRST [LSH G,1 ;[700] NO, IFF LO HALF=0!<377(OCT). CHANGE INDEX TO POINTER
;PICK UP DOUBLE WORD ANSWER
DMOVE A,POWERS(G)
JRST SCALE] ;SCALE RESULTS AND EXIT
;**; [700] ADD LABEL @ REDUCE+4L SWG 30-AUG-77
REDUC2: CAMGE B, TABLE(G) ;[700]IS ARGUMENT GREATER THAN ENTRY?
SOJA G, REDUCE ;NO, TRY NEXT LOWER ENTRY
SUB B, TABLE(G) ;YES, ALL DONE -REDUCE ARGUMENT
LSH G, 1 ;SAVE INDEX AS A POINTER
ASHC B, -8 ;MAKE ROOM FOR EXPONENT
IFE CPU-KA10,<MOVE A, B ;SET UP ARG. FOR NORMALIZING
ASH C, -8 ;MAKE ROOM FOR LOW ORDER EXP.
FSC A,200 ;SET EXP TO 200
FSC C,200-^D27 ;SET EXP 27 LOWER
FADL A,C > ;MAKE STANDARD NUMBER
IFE CPU-KI10,<TLO B,200000 ;INSERT EXPONENT
DFAD B,[EXP 0,0] ;NORMALIZE RESULT
DMOVE A,B > ;PUT RESULT AN A,B
;GET B4/F
DMOVE D,B4
FLDIV D,A
;SAVE B4/F
DMOVEM D,B4F
;GET F+B4/F
FLADD D,A
;GET H4/(F+B4/F)
DMOVEM D,FB4F
DMOVE D,H4
FLDIV D,FB4F
;GET D4*F
FLMUL A,D4
;GET (B4/F)-C4+D4*F+(H4/(F+B4/F))
FLADD D,A
FLADD D,MC4
FLADD D,B4F
;GET 1.0+A4/REST
DMOVE A,A4
FLDIV A,D
FLADD A,ONE
JUMPE G,SCALE ;MULTIPLY BY POWER OF TWO?
FLMUL A,POWERS(G)
SCALE: FSC A,.-. ;SCALE RESULTS
IFE CPU-KA10,<FSC B,@SCALE
JFCL ;SUPPRESS UNDERFLOW MSG. FROM LOW WORD.
FADL A,B > ;MAKE STANDARD NUMBER
DEXEND: MOVS X, XBLT ;PICK UP THE BLT POINTER
BLT X, X ;RESTORE THE ACCUMULATORS
GOODBY (1) ;EXIT
LIT
PRGEND
TITLE DLOG %4.(235) PDP-10/10I DOUBLE PRECISION LOGARITHM FUNCTION
SUBTTL H. P. WEISS/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1973,1977 BY DIGITAL EQUIPMENT CORPORATION
SEARCH FORPRM
ENTRY DLOG
EXTERN DLOG.
DLOG=DLOG.
PRGEND
TITLE DLOG10 %4.(235) PDP-10/10I DOUBLE PRECISION LOGARITHM FUNCTION
SUBTTL H. P. WEISS/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1973,1977 BY DIGITAL EQUIPMENT CORPORATION
ENTRY DLOG10
EXTERN DLG10.
DLOG10=DLG10.
PRGEND
TITLE DLOG. %4.(356) PDP-10/10I DOUBLE PRECISION LOGARITHM FUNCTION
SUBTTL D. TODD /DRT/HPW/MD 12-AUG-74
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1972,1977 BY DIGITAL EQUIPMENT CORPORATION
;THIS PROGRAM CALCULATES THE LOGARITHM OF A DOUBLE PRECISION
;ARGUMENT. THE ALGORITHM USED IS DESCRIBED ON PAGES 29-30 OF
;RALSTON AND WILF, "MATHEMATICAL METHODS FOR DIGITAL COMPUTERS".
;THE ARGUMENT X IS WRITTEN AS
; X = (2**N)*F WHERE 1/2 < F < 1
;THEN LOG(X) = (N*LOGE(2)) + LOG(F)
;F IS REDUCED BY FIXED POINT MULTIPLICATION BY NOT MORE THAN
;THREE CONSTANTS. THIS YIELDS
; 0 < T = A1*A2*A3*F - 1.0 < (2**-7)/5
;NOTE THAT NOT ALL THE A1,A2,A3 NEED BE INCLUDED IN THE PRODUCT.
;FINALLY,
; LOG(F) = LOG(1+T) - LOG(A1) - LOG(A2) - LOG(A3)
;LOG(1+T) IS CALCULATED AS A TAYLOR SERIES IN T.
;THE ROUTINE IS CALLED IN THE FOLLOWING MANNER:
; JSA Q, DLOG
; EXP ARG
;WHERE ARG IS THE ADDRESS OF THE HIGH ORDER PART OF THE DOUBLE
;PRECISION ARGUMENT. THE RESULT IS LEFT IN ACCUMULATOR A AND B.
;AVAILABLE TABLES FOR LOGARITHMS LIST 16 SIGNIFICANT DIGITS,
;ALL OF WHICH AGREE WITH VALUES PRODUCED BY THIS PROGRAM. AN
;ESTIMATE OF THE MAXIMUM THEORETICAL ERROR IS OBTAINED BY OB-
;SERVING THAT, AFTER REDUCTION, THE ARGUMENT IS LESS THAN
;(2**-7)/5 = .00625. HENCE, THE ERROR IN THE TAYLOR SERIES
;IS LESS THAN
; (.00625**9)/9 = 1.4573*10**-20/9 = .1619*10**-20
A= 0
B= 1
C= 2
D= 3
E= 4
F= 5
G= 6
Q= 16
P= 17
X= G ;HIGHEST AC TO SAVE
SEARCH FORPRM
;CONSTANTS AND TEMPORARY LOCATIONS AND STUFF
LOGLST: DOUBLE 200471174064,325425031470 ;0.61180 15411 05992 8976
DOUBLE 200402252251,151350376610 ;0.50455 60107 52395 2859
DOUBLE 177637144373,057714113734 ;0.40546 51081 08164 3810
DOUBLE 177506061360,207057302360 ;0.31845 37311 18534 6147
DOUBLE 176710776761,346515041520 ;0.22314 35513 14209 7553
DOUBLE 176537746034,051711723600 ;0.17185 02569 26659 2214
DOUBLE 175557032242,271265512760 ;0.08961 21586 89687 12374
DOUBLE 173770123303,236031377700 ;0.03077 16586 66753 68689
LIST: 354000000000 ;1.11011 BINARY
324000000000 ;1.10101 BINARY
300000000000 ;1.10000 BINARY
260000000000 ;1.01100 BINARY
240000000000 ;1.01000 BINARY
230000000000 ;1.00110 BINARY
214000000000 ;1.00011 BINARY
204000000000 ;1.00001 BINARY
DLOGE2: DOUBLE 200542710277,276434757153 ;0.69314 71805 59945 30941 72321
A9: DOUBLE 175707070707,034343434344 ;1/9
A8: DOUBLE 601400000000,0 ;-1/8
A7: DOUBLE 176444444444,222222222222 ;1/7
A6: DOUBLE 601252525252,652525252526 ;-1/6
A5: DOUBLE 176631463146,146314631463 ;1/5
A4: DOUBLE 600400000000,0 ;-1/4
A3: DOUBLE 177525252525,125252525253 ;1/3
A2: DOUBLE 577400000000,0 ;-1/2
ONE:
A1: 201400000000
DZERO: 000000000000
IFE CPU-KI10,<0 >
SUMSAV: BLOCK 2 ;STORAGE FOR PARTIAL ANSWER
XBLT: XWD C, ACSAVE
ACSAVE: BLOCK X-C+1
;DOUBLE PRECISION COMMON LOGARITHM FUNCTION
;THE ROUTINE CALCULATES THE LOGARITHM, BASE 10, OF A DOUBLE
;PRECISION ARGUMENT. THE ALGORITHM USED IS
; DLOG10(X) = DLOG(X)*LOG10(E)
;THE CALLING SEQUENCE FOR THE ROUTINE IS THE FOLLOWING:
; JSA Q, DLOG10
; EXP ARG
;THE DOUBLE PRECISION ANSWER IS RETURNED IN AC A AND B.
HELLO (DLG10.,DLOG10) ;[235] ENTRY TO DLOG10 ROUTINE
MOVEI B,@(Q) ;GET ADDRESS OF ARGUMENT
SKIPN B,(B) ;IF ARG = 0,
JRST ZERO ;THEN GO TO ZERO.
PUSHJ P,DLOG. ;CALCULATE LOG(ARG)
IFE CPU-KA10,<MOVEM C,ACSAVE > ;SAVE ACCUMULATOR C
FLMUL A,LOG10D
IFE CPU-KA10,<MOVE C,ACSAVE > ;RESTORE C
GOODBY (1) ;EXIT
LOG10D: DOUBLE 177674557305,111562416145 ;.43429448190325182765
HELLO (DLOG,.) ;[235] ENTRY TO DLOG ROUTINE
MOVE 0,XBLT ;SAVE AC'S
BLT 0,ACSAVE+X-C ;...
DMOVE A,@(Q) ;PICK UP ARGUMENT
JUMPG A,DLOGP ;ARGUMENT IS GREATER THAN 0
JUMPE A, ZERO ;CHECK FOR ZERO ARGUMENT
ERROR (LIB,12,2,[ASCIZ /Attempt to take DLOG of Negative Arg/])
DMOVN A,@(Q) ;GET ABSF(X)
DLOGP: CAMN A,ONE ;X PRECISELY 1.0?
JUMPE B,[SETZB A,B
JRST DLOG5] ;YES, RETURN ZERO
DLOG1: LDB D,[POINT 8,A,8] ;NO,PICK UP EXPONENT FROM HI ORDER
SUBI D, 200 ;GET EXPONENT EXCESS 200
FSC D,233 ;MAKE FLOATING POINT NUMBER
MOVEI E,0 ;SET UP LO HALF
;CALCULATE N*LOGE(2)
FLMUL D,DLOGE2
DMOVEM D,SUMSAV
IFE CPU-KA10,<LSH B,8 > ;GET RID OF LOW ORDER EXP.
TLZ A,777000 ;MASK OUT EXPONENT
ASHC A,8 ;NORMALIZE FRACTION TO BIT 1
SETZB D,E ;INITIALIZE REDUCTION CONSTANT TO 0
DLOG2: LDB G,[POINT 3,A,4] ;GET HI 3 BITS BELOW 1/2
MUL B, LIST(G) ;FIXED POINT MULTIPLY FOR LO HALF
MOVE C,B ;SAVE HI HALF OF LO PRODUCT
;(LO HALF IS ALL 0'S, THROW IT AWAY)
MUL A, LIST(G) ;MULTIPLY HI ORDER, TOO
TLO B,(1B0) ;SET BIT 0, TO AVOID OVERFLOW
ADD B,C ;COMBINE RESULTS OF MULTIPLY
TLZN B,(1B0) ;CLEAR BIT 0, WAS THERE OVERFLOW?
ADDI A,1 ;YES, PROPOGATE CARRY
ASH G, 1 ;TURN BITS INTO D.P. POINTER
FLADD D,LOGLST(G)
TLZE A, 200000 ;IS THE PRODUCT .GE. 1.0?
JRST DLOG3 ;YES
ASHC A, 1 ;NO, GET RID OF EXTRANEOUS ZERO
JRST DLOG2 ;TRY ANOTHER MULTIPLICATION
DLOG3: ASHC A, -7 ;MAKE ROOM FOR THE EXPONENT
IFE CPU-KA10,<FSC A,200 ;SET EXPONENT TO 200
ASH B,-8 ;MAKE ROOM FOR LO EXPONENT
FSC B,200-^D27 ;INSERT LO EXPONENT
FADL A,B > ;MAKE NORMAL DOUB, PRECISION NUMBER
IFE CPU-KI10,<TLC A,200000 ;INSERT EXPONENT
DFAD A,DZERO > ;NORMALIZE
DFN D,E ;NEGATE LOG OF MAGIC NUMBERS
;AND ADD INTO FINAL SUMMATION
FLADD D,SUMSAV
DMOVEM D,SUMSAV
;PICK UP CONSTANT TO START
DMOVE D,A9
MOVEI G,A8 ;INIT TABLE POINTER AT A8
DLOG4: ;MULTIPLY ACCUMULATED SUM BY X
FLMUL D,A
;ADD NEXT CONSTANT INTO PARTIAL SUM
FLADD D,0(G)
ADDI G, 2 ;UPDATE POINTER TO NEXT CONSTANT
CAIG G, A1 ;ARE WE DONE YET?
JRST DLOG4 ;NO, LOOP BACK FOR MORE TAYLOR SERIES
;YES, ONE LAST MULTIPLICATION
FLMUL A,D
;AND ADD SERIES SUM INTO FINAL ANSWER
FLADD A,SUMSAV
DLOG5: MOVS X, XBLT ;PICK UP BLT POINTER
BLT X,X ;RESTORE ACCUMULATORS
GOODBY (1) ;EXIT
ZERO: MOVSI A,400000 ;SET UP LARGE NEG. NUM.
IFE CPU-KA10,<FADL A,A > ;CAUSE NEGATIVE OVERFLOW
IFE CPU-KI10,<DFAD A,A >
GOODBY (1) ;TRAP ROUTINE FIXED IT UP AND
;PRINTED MESSAGE
PRGEND
TITLE DSQRT %4.(235) DOUBLE PRECISION SQUARE ROOT
SUBTTL H. P. WEISS/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1973,1977 BY DIGITAL EQUIPMENT CORPORATION
SEARCH FORPRM
ENTRY DSQRT
EXTERN DSQRT.
DSQRT=DSQRT.
PRGEND
TITLE DSQRT. %4.(356) DOUBLE PRECISION SQUARE ROOT
SUBTTL D. TODD /DRT/HPW/MD 12-AUG-74
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1972,1977 BY DIGITAL EQUIPMENT CORPORATION
;FROM V.022 31-DEC-69
;FROM V.020 5 MAY, 1969 /TWE
;FROM V.005 2-MAR-67
;DOUBLE PRECISION SQUARE ROOT FUNCTION
;THIS ROUTINE CALCULATES THE SQUARE ROOT OF A DOUBLE PRECISION
;ARGUMENT BY DOING A LINEAR SINGLE PRECISION APPROXIMATION ON
;THE HIGH ORDER WORD, THEN TWO DOUBLE PRECISION ITERATIONS OF
;NEWTONS METHOD. THIS SHOULD GENERATE A RESULT ACCURATE TO
;20 DECIMAL SIGNIFICANT DIGITS. THE ALGORITHM IS AS FOLLOWS
;X = (2**(2N))*F, WHERE 1/2 < F < 1
;HENCE SQRT(X) = 2**N * SQRT(F)
;THE LINEAR APPROXIMATION IS OF THE FORM
;SQRT(F) = A2 - B2/(C2+F-D2/(E2+F))
;WHERE THE CONSTANTS A2,B2,C2,D2, AND E2 HAVE THE FOLLOWING
;VALUES
;CONSTANT VALUE WHEN 0.25<F<0.50 VALUE WHEN 0.50<F<1.0
;A2 (5/14)*SQRT(70) (5/7)*SQRT(35)
;B2 (50/49)*SQRT(70) (200/49)*SQRT(35)
;C2 47/14 47/7
;D2 4/49 16/49
;E2 3/14 3/7
A= 0
B= 1
C= 2
D= 3
E= 4
F= 5
Q= 16
P= 17
X= F ;HIGHEST AC SAVED
SEARCH FORPRM
;CONSTANTS AND TEMPORARY LOCATIONS AND STUFF
A2: 202576362203 ;2.98807152
203416346045 ;4.225771271
B2: 204421143713 ;8.537347194
205602266310 ;24.14726441
C2: 202655555556 ;3.357142857
203655555556 ;6.7142857143
D2: 175516274052 ;0.0816326531
177516274052 ;0.326530612
E2: 176666666667 ;0.2142857143
177666666667 ;0.4285714286
XBLT: XWD C,ACSAVE
ACSAVE: BLOCK X-C+1
N: BLOCK 2 ;STORAGE FOR ARGUMENT
HELLO (DSQRT,.) ;[235] ENTRY TO DSQRT ROUTINE
MOVE 0, XBLT ;SET UP BLT POINTER
BLT 0, ACSAVE+X-C ;SAVE SOME ACCUMULATORS
DMOVE A,@(Q) ;GET D.P. ARG
JUMPG A, DSQRTP ;ARGUMENT IS GREATER THAN 0
JUMPE A, DSQRT4 ;ARGUMENT OF ZERO?
;**; [554] CHANGE @ N + 7L IN DSQRT. CLRH 11-JUN-76
ERROR (LIB,13,2,[ASCIZ /Attempt to take DSQRT of Negative Arg/])
DMOVN A,@(Q) ;GET ARGS AGAIN (ABSF)
DSQRTP: MOVE F, A ;GET SPARE COPY OF HIGH ORDER
LSH F, -33 ;GET RID OF FRACTION BITS
SUBI F, 201 ;GET RID OF THE BASE 200 PART OF
;EXPONENT. EXTRA 1 IS A FUDGE.
ROT F, -1 ;CUT EXPONENT IN HALF, SAVE EXTRA
;BIT FOR LATER USE AS INDEX REG.
HRRM F, DSQRT1 ;SAVE REDUCED EXPONENT FOR SCALING
LSH F, -43 ;BRING BIT BACK - IF 0, THEN
;1/4<F<1/2,OTHERWISE 1/2<F<1.
TLZ A, 777000 ;WIPE OUT EXPONENT BITS IN ARG.
FSC A, 177(F) ;RESET IT TO EITHER 177 OR 200
IFE CPU-KA10,<TLZ B,777000 ;WIPE OUT EXP BITS IN LOW ARG
FSC B,177-^D27(F) ;SET IT TO 27 LESS THAN HI PART
FADL A,B > ;UNNORMALIZE LOW PART
MOVE D, A ;PICK UP ANOTHER COPY OF NEW FRAC.
FADR D, E2(F) ;FORM E2+F
MOVN C, D2(F) ;PICK UP -D2
FDVR C, D ;CALCULATE -D2/(E2+F)
FADR C, C2(F) ;GET C2-D2/(E2+F)
FADR C, A ;CALCULATE F+C2-D2/(E2+F)
MOVN D, B2(F) ;PICK UP -B2
FDVR D, C ;GET -B2/(F+C2-D2/(E2+F))
FADR D, A2(F) ;GET FINAL FIRST APPROXIMATION
MOVEI E,0 ;LOW HALF OF 1ST APPROX. IS 0
DMOVEM A,N ;SAVE DSQRT ARGUMENT
;GET N/X0
FLDIV A,D
;X0+N/X0
FLADD D,A
FSC D,-1 ;X1=.5*(X0+N/X0)
IFE CPU-KA10,<FSC E,-1 ;...
FADL D,E > ;UNNORMALIZE LOW WORD
DMOVE A,N ;GET N BACK
;N/X1
FLDIV A,D
;X1+N/X1
FLADD A,D
DSQRT1: FSC A,.-. ;SCALE RESULTS FOR ANSWER
IFE CPU-KA10,<FSC B,@DSQRT1 ;...
FADL A,B > ;UNNORMALIZE ANSWER
DSQRT3: MOVS X, XBLT ;SET UP THE BLT POINTER
BLT X, X ;RESTORE THE ACCUMULATORS
DSQRT4: GOODBY (1) ;EXIT
PRGEND
TITLE DSIN %4.(235) DOUBLE PRECISION SIN AND COSINE ROUTINE
SUBTTL H. P. WEISS/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1973,1977 BY DIGITAL EQUIPMENT CORPORATION
SEARCH FORPRM
ENTRY DSIN
EXTERN DSIN.
DSIN=DSIN.
PRGEND
TITLE DCOS %4.(235) DOUBLE PRECISION SIN AND COSINE ROUTINE
SUBTTL H. P. WEISS/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1973,1977 BY DIGITAL EQUIPMENT CORPORATION
ENTRY DCOS
EXTERN DCOS.
DCOS=DCOS.
PRGEND
TITLE DSIN. %4.(356) DOUBLE PRECISION SIN AND COSINE ROUTINE
SUBTTL D. TODD /DRT/HPW/MD 12-AUG-74
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1972,1977 BY DIGITAL EQUIPMENT CORPORATION
;DOUBLE PRECISION SINE AND COSINE FUNCTION
;THIS PROGRAM CALCULATES THE SINE AND COSINE OF A DOUBLE
;PRECISION ARGUMENT. FOR A DESCRIPTION OF THE ALGORITHM,
;SEE THE SCIENCE LIBRARY AND FORTRAN UTILITY SUBPROGRAMS
;MANUAL.
;THE CALLING SEQUENCE FOR THE PROGRAM IS AS FOLLOWS:
; JSA Q, DSIN
; EXP ARG
;WHERE ARG IS THE ADDRESS OF THE HIGH ORDER PORTION OF THE
;ARGUMENT. THE DOUBLE PRECISION ANSWER IS LEFT IN ACCUMULATORS
;A AND B.
SEARCH FORPRM
A= 0
B= 1
C= 2
D= 3
E= 4
F= 5
G= 6
Q= 16
P= 17
X=G ;SAVE AC'S C THRU G
;CONSTANTS AND TEMPORARY LOCATIONS AND STUFF
C17: DOUBLE 120625130734,014126512326 ;1/17!=.28114572543455207632E-14
DOUBLE 124656376371,314734037043 ;1/16!=.47794773323873852974E-13
C15: DOUBLE 647121401406,463043740735 ;-1/15!=-.76471637318198164759E-12
DOUBLE 643154321325,717701542677 ;-1/14! =-.11470745597729724714E-10
C13: DOUBLE 140541110604,352066411370 ;1/13!=.16059043836821614599E-9
DOUBLE 144436733073,376154227552 ;1/12!=.20876756987868098979E-8
C11: DOUBLE 630121467246,402535434340 ;-1/11!=-.25052108385441718775E-7
DOUBLE 624330066022,441660243433 ;-1/10!=-.27557319223985890653E-6
C9: DOUBLE 156561674351,125543463437 ;1/9!=.27557319223985890653E-5
DOUBLE 161640064006,200320032003 ;1/8!=.24801587301587301587E-4
C7: DOUBLE 613137713771,577457745775 ;-1/7!=-.19841269841269841270E-3
DOUBLE 610223722372,517511751175 ;-1/6!=-.1388888888888888889E-2
C5: DOUBLE 172421042104,104210421042 ;1/5!=.00833333333333333333333
DOUBLE 174525252525,125252525253 ;1/4!=.041666666666666666667
C3: DOUBLE 601252525252,652525252526 ;-1/3!=-0.16666666666666666667
C2: 577400000000 ;-1/2!=-0.50000000000000000000
000000000000
PIOTLO=021026430215 ;LOW HALF OF PI/2 FOR PDP-6/KI10
ADDTBL: DOUBLE 576155700452,-PIOTLO ;-PI/2
DOUBLE 575155700452,-PIOTLO ;-PI
574322320340 ;-3*PI/2
IFE CPU-KA10,<150146336134>
IFN CPU-KA10,<463157055627>
DOUBLE 574155700452,-PIOTLO ;-2*PI
ONE: 1.0
0
XBLT: XWD C, ACSAVE
PIOT: DOUBLE 201622077325,PIOTLO ;PI/2
FLAG: 0 ;NEGATE ANSWER IF (FLAG) IS MINUS
COSFLG: BLOCK 1 ;2 FOR COSINE SERIES, 0 FOR SINE
ACSAVE: BLOCK X-C+1 ;SAVE ROOM FOR STORING AC'S
ARGSAV: BLOCK 2 ;TEMP FOR X
HELLO (DCOS,.) ;[235] ENTRY TO COSINE ROUTINE
MOVE 0,XBLT
BLT 0,ACSAVE+X-C ;SAVE ACCUMULATORS
DMOVE A,@(Q) ;GET DOUBLE WORD ARGUMENT
FLADD A,PIOT ;COS(X)=SIN(PI/2+X), LEAVE RESULT IN A
JRST DSIN0 ;ENTER SINE ROUTINE
HELLO (DSIN,.) ;[235] ENTRY TO SINE ROUTINE
MOVE 0,XBLT
BLT 0,ACSAVE+X-C ;SAVE ACCUMULATORS
DMOVE A,@(Q) ;GET DOUBLE WORD ARGUMENT
DSIN0: JUMPE A, DSIN4C ;ARGUMENT OF ZERO?
SETZB G, FLAG ;SET FLAG FOR POSITIVE ARGUMENT
JUMPGE A, DSIN1 ;IS ARGUMENT POSITIVE?
SETOM FLAG ;NO, CHANGE FLAG
DFN A, B ;NEGATE THE ARGUMENT
DSIN1: DMOVEM A,ARGSAV
;CALCULATE X/(PI/2),LEAVE THE RESULTS IN A,A+1
FLDIV A,PIOT
CAML A,[1.0] ;IS X .L. PI/4?
JRST REDUCE ;NO, REDUCE IT
CAML A,[0.5] ;IS X .GE. PI/4
MOVEI G,1 ;YES, 2ND OCTANT
DSIN1A: DMOVE A,ARGSAV
DSIN2: TRNE G,4 ;QUADRANTS 3 OR 4?
SETCMM FLAG ;YES, SINE IS NEGATIVE
JUMPE G,DSIN3 ;IS X .L. PI/4
TRNE G,1 ;NO,GET INDEX INTO QUADRANT TABLE
ADDI G,1 ;...
;MAKE: -PI/4 .GE. X .LE. +PI/4
FLADD A,ADDTBL-2(G)
SKIPGE A
DFN A,B ;TAKE ABSOLUTE VALUE
DMOVEM A,ARGSAV
DSIN3: TRZ G,777775 ;LEAVE ONLY OCTANT BIT
HRRZM G,COSFLG ;0 FOR SINE SERIES, 2 FOR COSINE
CAMG A,[XWD 147471,421605] ;IS X .L. SQRT(6)*2**-27?
JRST SMALL ;YES, THEN SIN(X)=X
;CALCULATE X**2, AND LEAVE IN A,A+1
FLMUL A,ARGSAV
;INITIALIZE PARTIAL SUM
DMOVE D,C17(G)
MOVEI G,C15(G) ;TURN OCTANT POINTER INTO TABLE ADR
DSIN3C: ;MULTIPLY PARTIAL SUM BY X**2
FLMUL D,A
;ADD NEXT CONSTANT TO PARTIAL SUM
FLADD D,0(G)
ADDI G,4 ;MOVE POINTER TO NEXT CONSTANT
CAIG G,C2 ;DONE?
JRST DSIN3C ;NO, LOOP BACK FOR MORE OF SERIES
;YES, ONE MORE MULTIPLY
FLMUL A,D
;ADD 1.0 INTO SUM
FLADD A,ONE
SKIPE COSFLG ;IS THIS COSINE SERIES?
JRST DSIN4A ;YES
;NO, MULTIPLY BY X, THIS IS SIN
FLMUL A,ARGSAV
DSIN4A: SKIPE FLAG ;NEGATE RESULT?
DFN A,B ;YES
DSIN4C: MOVS X,XBLT ;SET UP BLT POINTER
BLT X,X ;RESTORE AC'S
GOODBY (1) ;EXIT
SMALL: JUMPE G,DSIN4A ;CALCULATING COSINE?
DMOVE A,ONE ;YES, COS(X)=1.0
JRST DSIN4A
REDUCE: MOVE D,A ;SAVE QUADRANT NUMBER
LDB G,[POINT 8,A,8] ;GET EXPONENT
IFE CPU-KA10,<
LSH B,9> ;WIPE OUT LO EXPONENT
IFN CPU-KA10,<
LSH B,1> ;REMOVE 1 BIT FOR KI10
TLZ A,777000 ;DITTO HI EXPONENT
LSHC A,-202(G) ;MAKE ARGUMENT MODULO 2 PI
LDB G,[POINT 3,A,11] ;GET QUADRANT AND OCTANT BITS
CAMGE D,[4.0] ;IS NON-REDUCED ARG OK?
JRST DSIN1A ;YES, SAVE THE DFMP INACCURICIES
TLZ A,777000 ;MAKE WAY FOR EXPONENT
IFE CPU-KA10,<
FSC A,202 ;PUT EXP IN HI WORD
LSH B,-9 ;MAKE WAY FOR LOW EXPONENT
FSC B,202-^D27 ;PUT IN LO EXPONENT
FADL A,B> ;UNNORMALIZE LO WORD
IFN CPU-KA10,<
TLO A,202000 ;PUT EXP IN HI WORD, NO NORMALIZE
LSH B,-1 ;PUT LO FRACTION IN PLACE
DFAD A,[EXP 0, 0]> ;NORMALIZE
FLMUL A,PIOT ;CHANGE MAKE TO RADIANS (MOD 2 PI)
DMOVEM A,ARGSAV
JRST DSIN2 ;GO CHANGE ARGUMENT TO 1ST OCTANT
LIT
PRGEND
TITLE DATAN2 %4.(235) TWO ARGUMENT DOUBLE PRECISION ARC TANGENT
SUBTTL H. P. WEISS/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1973,1977 BY DIGITAL EQUIPMENT CORPORATION
SEARCH FORPRM
ENTRY DATAN2
EXTERN DATN2.
DATAN2=DATN2.
PRGEND
TITLE DATN2. %4.(356) TWO ARGUMENT DOUBLE PRECISION ARC TANGENT
SUBTTL D. TODD /DRT/MD 12-AUG-74
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1972,1977 BY DIGITAL EQUIPMENT CORPORATION
;THIS ROUTINE CALCULATES THE ARCTANGENT OF A/B
;IF ARGUMENT IS IN 2ND QUADRANT, DATAN2(A,B)=DATAN(A/B) + PI
;IF ARGUMENT IS IN 3RD QUADRANT, DATAN2(A,B)=DATAN(A/B)-PI
;IF ARGUMENT IS IN 1ST OR 4TH QUADRANT, DATAN2(A,B)=DATAN(A/B)
;IF QUOTIENT A/B OVER OR UNDERFLOWS, RETURN AN ANGLE
;ON A CO-ORDINATE AXIS
;THE ROUTINE IS CALLED IN THE FOLLOWING MANNER:
; JSA Q, DATAN2
; EXP A
; EXP B
;DATAN(A/B) IS RETURNED AS A DOUBLE PRECISION ANSWER IN ACCUMULATORS
;A AND B.
SEARCH FORPRM
A= 0
B= 1
C= 2
D= 3
Q= 16
P= 17
;CONSTANTS AND TEMPORARY LOCATIONS AND STUFF
DARG: BLOCK 2 ;ARGUMENT FOR DATAN
IFE CPU-KA10,<
CSAVE: BLOCK 1 >
PIOT: DOUBLE 201622077325,021026430215 ;+PI/2
PI.: DOUBLE 202622077325,021026430215 ;DEC=3.14159265358979323846
HELLO (DATN2.,DATAN2) ;[235] ENTRY TO DATAN2 ROUTINE
DMOVE A,@(Q) ;PICK UP 1ST ARG
IFE CPU-KA10,<MOVEM C,CSAVE ;SAVE AC C
MOVEM B,DARG ;SAVE LOW PART OF ARG.
MOVEI B,@1(Q) ;PICK UP
MOVE C,(B) ;THE
MOVE B,1(B) ;NEXT ARG.
EXCH B,DARG ;RESTORE LOW OF ONE ARG AND
MOVEM C,DARG+1 ;COPY THE OTHER ARG.
FDVL A,C ;CALCULATE A/B
JFCL AXIS ;UNDER/OVER FLOW
MOVN C,A ;...
FMPR C,DARG ;...
JFCL
UFA B,C
FDVR C,DARG+1
JFCL
FADL A,C
JFCL AXIS > ;UNDER/OVER FLOW
IFE CPU-KI10,<DFDV A,@1(Q) ;CALCULATE A/B
JFCL AXIS > ;TRANSFER FOR OVER/UNDER FLOW
DMOVEM A,DARG ;STORE ARG FOR DATAN
FUNCT DATAN.,<DARG> ;CALCULATE DATAN(A/B)
IFE CPU-KA10,<MOVE C, CSAVE > ;RESTORE ACCUMULATOR C
SKIPL @1(Q) ;WAS 2ND ARGUMENT POSITIVE?
GOODBY (2) ;YES, 1ST OR 4TH QUAD, EXIT
;NO, 2ND OR 3RD
SKIPGE @(Q) ;2ND QUADRANT?
DFN A,B ;NO, 3RD. SUBTRACT PI
FLADD A,PI.
IFE CPU-KA10,<MOVE C,CSAVE > ;RESTORE AC C.
SKIPGE @(Q) ;3RD OR 4TH QUADRANTS?
DFN A,B ;YES, NEGATE FINAL ANSWER
GOODBY (2) ;EXIT
AXIS:
IFE CPU-KA10,<MOVE C,CSAVE > ;RESTORE AC C
JUMPN A,OVER ;GO TO OVER IF OVERFLOW.
SKIPL @1(Q) ;ANS UNDERFLOWS IF Y/X UNDERFLOWS
JRST UNDMSG ;AND IF X >= 0.
DMOVE A,PI. ;O'E, ANS = +-PI, SO
JRST SETSGN ;GO TO SET SIGN.
OVER: SKIPN @(Q) ;WAS Y =0 ?
FDVR 0,@(Q) ;YES, THIS IS 0/0 , FORCE A DIV CHK MSG.
DMOVE A,PIOT ;ANS = +-PI/2.
SETSGN: SKIPGE @(Q) ;ANS > 0 IF Y > 0.
DFN A,B ;ANS < 0 IF Y > 0.
GOODBY (2) ;EXIT.
UNDMSG: ERROR (APR,7,1,.+1) ;RETURN UNDERFLOW
SETZB 0,1 ;AND ANS = 0
GOODBY (2) ;EXIT.
PRGEND
TITLE DATAN %4.(235) SINGLE ARGUMENT DOUBLE PRECISION ARC TANGENT
SUBTTL H. P. WEISS/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1973,1977 BY DIGITAL EQUIPMENT CORPORATION
SEARCH FORPRM
ENTRY DATAN
EXTERN DATAN.
DATAN=DATAN.
PRGEND
TITLE DATAN. %4C.(513) SINGLE ARGUMENT DOUBLE PRECISION ARC TANGENT
SUBTTL D. TODD /DRT/HPW/MD/JNG 15-DEC-75
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1972,1977 BY DIGITAL EQUIPMENT CORPORATION
;THIS ROUTINE CALCULATES THE ACTANGENT OF A DOUBLE PRECISION
;ARGUMENT ACCORDING TO THE ALGORITHM
;DATAN(X) = LAMBDA*X/(Z+B0+A1/(Z+B1+A2/(Z+B2+A3/(Z+B3))))
;FOR X>1.0, THE IDENTITY
; ATAN(X) = PI/2 - ATAN(1/X)
;IS USED. FOR 0.236<X<1.0, THE IDENTITY
; ATAN(X) = ATAN(1/2) + ATAN(2X-1/X+2)
;IS USED.
;FOR X<SQRT(3)*2**-27, ATAN(X)=X IS USED
;THE ROUTINE HAS THE FOLLOWING CALLING SEQUENCE
; JSA Q, DATAN
; EXP ARG
;WHERE ARG IS THE ADDRESS OF THE HIGH ORDER PART OF THE DOUBLE
;PRECISION WORD. THE ANSWER IS RETURNED IN ACCUMULATORS A AND B.
SEARCH FORPRM
A= 0
B= 1
C= 2
D= 3
E= 4
F= 5
G= 6 ;FLAG REGISTER
;BIT35=1, ADD ATAN(1/2) TO ANSWER
;BIT34=1, ADD 2*ATAN(1/2) TO ANSWER ;[513]
;BIT17=1,ADD -PI/2 TO ANSWER ;[513]
;BIT0=1, NEGATE FINAL ANSWER
H= 7 ;LOOP POINTER
Q= 16
P= 17
X= H ;NUMBER OF ACCUMULATORS TO SAVE
;REVISION HISTORY
;EDIT SPR DESCRIPTION
;---- --- -----------
;513 15636 WHEN (5*SQRT(5)-2)/11 < ABS(X) < (5*SQRT(5)+2)/11,
; THE DATAN(X) RETURNED CAN BE WRONG AFTER THE 13TH PLACE.
; THIS IS BECAUSE (2*X-1)/(X+2) > SQRT(5)-2 IN THIS RANGE.
;CONSTANTS AND TEMPORARY LOCATIONS AND STUFF
LAMBDA: DOUBLE 204613772770,017027645561 ;12.37469 38775 51020 40816
B0: DOUBLE 205644272446,121335250615 ;26.27277 52490 26980 67155
A1: DOUBLE 570276502107,437176661671 ;-80.34270 56102 16599 70467
B1: DOUBLE 203627237361,165414142742 ;6.36424 16870 04411 34492
A2: DOUBLE 576316772502,512470127251 ;-1.19144 72238 50426 48905
B2: DOUBLE 202415301602,015271031674 ;2.10451 89515 40978 95180
A3: DOUBLE 602277106546,717167531241 ;-0.07833 54278 56532 11777
B3: DOUBLE 201502125320,370207664057 ;1.25846 41124 27629 031727
ATANH: DOUBLE 177732614701,130335517321 ;ATAN(1/2)
MONE: EXP -1.0,0
TWO: EXP +2.0,0
MPIOT: DOUBLE 576155700452,756751347563 ;-PI/2
SMALL: XWD 146673,317272 ;SQRT(3)*2**-27
XBLT: XWD C, ACSAVE
ACSAVE: BLOCK X-C+1
DX: BLOCK 2 ;HOLDS X MOST OF THE TIME
X2: BLOCK 2 ;HOLDS X**2
HELLO (DATAN,.) ;[235] ENTRY TO DATAN ROUTINE
MOVE 0,XBLT ;SAVE ACCUMULATORS
BLT 0,ACSAVE+X-C ;...
DMOVE A,@(Q) ;GET ARGUMENT
JUMPE A,DATAN6 ;ARG .E. 0?
HLLZ G, A ;LH(G)=SGN(A), RH(G) = 0
;**;[513] INSERT @ DATAN.+7L JNG 15-DEC-75
TLZ G,377777 ;[513] ZAP ALL BUT SIGN FOR FLAGS
SKIPGE A ;IS THE ARGUMENT POSITIVE?
DFN A, B ;NO,NEGATE IT
MOVSI D, (1.0) ;GET DOUBLE PRECISION 1.0
MOVEI E,0 ;0 LO PART
CAMN A,D ;IS HIGH ORDER EQUAL TO 1.0?
SKIPE B ;YES, IS LOW ORDER ZERO?
CAMGE A,D ;NO, IS ARG>1.0?
JRST DATAN0 ;NO
TLC G,(1B0) ;COMPLEMENT FINAL SIGN BIT, GET 1/X
;**;[513] CHANGE @ DATAN.+17L JNG 15-DEC-75
TLO G,1 ;[513] ADD -PI/2 TO FINAL ANSWER
FLDIV D,A
DMOVEM D,A
DATAN0: DMOVEM A,DX
CAMGE A,[0.236] ;IS ARG .GE. (SQRT(5)-2)?
JRST DATAN1 ;NO, PROCEED WITH ALGORITHM
;CALCULATE X+2
FLADD A,TWO
EXCH A,DX ;GET X, SAVE X+2
EXCH B,DX+1 ;...
FSC A,1 ;CALCULATE 2X
IFE CPU-KA10,<FSC B,1 ;...
FADL A,B > ;...
;CALCULATE 2X-1
FLADD A,MONE
;(2X-1)/(X+2) WITH RESULTS IN A,B
FLDIV A,DX
;**;[513] REPLACE @ DATAN0+14L JNG 15-DEC-75
AOJA G,DATAN0 ;[513] TRY AGAIN IN CASE STILL TOO BIG
DATAN1: MOVM D,A ;GET MOD(X)
CAMGE D,SMALL ;CAN ATAN(X)=X?
JRST DATAN3 ;YES
;CALCULATE X**2
FLMUL A,DX
DMOVEM A,X2
;INIT CONTINUED FRACTION COMP. WITH B3
DMOVE A,B3
MOVEI H,B3 ;INIT POINTER TO NUMBER TABLE
JRST DAT2 ;DIVE INTO LOOP
DAT1: ;ADD B1
FLADD A,0(H)
DAT2: ;ADD X**2
FLADD A,X2
;GET A3 (OR A1)
DMOVE D,-2(H)
FLDIV D,A
;ADD B2 (OR B0)
FLADD D,-4(H)
;ADD X**2
FLADD D,X2
;GET A2 (OR LAMBDA)
DMOVE A,-6(H)
FLDIV A,D
SUBI H,8 ;DECREMENT TABLE POINTER
CAILE H,LAMBDA ;FINISHED?
JRST DAT1 ;NO, DO IT LAST TIME
;MULTIPLY BY X
FLMUL A,DX
DATAN3:
;**;[513] CHANGE @ DATAN3+1L JNG 15-DEC-75
IFE CPU-KA10,<TRNN G,-1 ;[513] ADD ATAN(1/2)?
JRST DATAN4 ;[513] NO
FLADD A,ATANH
SOJA G,DATAN3 ;[513] TRY AGAIN IN CASE NEED TO ADD
;[513] ATAN(1/2) TWICE
DATAN4: TLNN G,1 ;[513] ADD -PI/2?
JRST DATAN5 ;NO
FLADD A,MPIOT
>
DATAN5:
IFE CPU-KI10,<TRNN G,-1 ;[513] NEED TO ADD ATAN(1/2)?
JRST DATAN7 ;[513] NO, PROCEED
DFAD A,ATANH
SOJA G,DATAN5 ;[513] MAKE SURE ALL DONE
DATAN7: TLNE G,1 ;[513] NEED TO ADD -PI/2?
DFAD A,MPIOT >
SKIPGE G ;NEGATE RESULT?
DFN A,B ;YES
DATAN6: MOVS X,XBLT ;RESTORE ACCUMULATORS
BLT X,X ;...
GOODBY (1) ;EXIT
PRGEND
TITLE DFLOAT %4.(235) INTEGER TO DOUBLE PRECISION CONVERSION
SUBTTL H. P. WEISS/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1973,1977 BY DIGITAL EQUIPMENT CORPORATION
SEARCH FORPRM
ENTRY DFLOAT
EXTERN DFLOT.
DFLOAT=DFLOT.
PRGEND
TITLE DFLOT. %4.(235) INTEGER TO DOUBLE PRECISION CONVERSION
SUBTTL ED YOURDON/KK/VAA/TWE/DMN
SUBTTL D. TODD /DRT/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1972,1977 BY DIGITAL EQUIPMENT CORPORATION
;36 BIT FLOAT FUNCTION
;CONVERTS A SIGNED FIXED POINT INTEGER TO FLOATING POINT
;BY BREAKING THE INTEGER INTO HIGH ORDER AND LOW ORDER
;FRACTIONS, CALCULATING AN EXPONENT, THEN ADDING THE TWO
;TOGETHER. DOUBLE PRECISION CONVERSION.
;THE ROUTINE IS CALLED AS FOLLOWS:
; JSA Q,D FLOAT
; EXP ARG
;THE ANSWER IS RETURNED IN ACCUMULATORS A AND B
SEARCH FORPRM
A= 0
B= 1
Q= 16
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 %4.(235) DOUBLE PRECISION TO INTEGER CONVERSION
SUBTTL H. P. WEISS/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1973,1977 BY DIGITAL EQUIPMENT CORPORATION
SEARCH FORPRM
ENTRY IDINT
EXTERN IDINT.
IDINT=IDINT.
PRGEND
TITLE DFIX %4.(235) DOUBLE PRECISION TO INTEGER CONVERSION
SUBTTL H. P. WEISS/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1973,1977 BY DIGITAL EQUIPMENT CORPORATION
SEARCH FORPRM
ENTRY DFIX
EXTERN DFIX.
DFIX=DFIX.
PRGEND
TITLE IDFIX %4.(235) DOUBLE PRECISION TO INTEGER CONVERSION
SUBTTL H. P. WEISS/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1973,1977 BY DIGITAL EQUIPMENT CORPORATION
SEARCH FORPRM
ENTRY IDFIX
EXTERN IDFIX.
IDFIX=IDFIX.
PRGEND
TITLE IDINT. %4.(235) DOUBLE PRECISION TO INTEGER CONVERSION
SUBTTL TOM EGGERS/DMN
SUBTTL D. TODD/DRT/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1972,1977 BY DIGITAL EQUIPMENT CORPORATION
Q=16
SEARCH FORPRM
HELLO (IDFIX,.) ;[235] ENTRY TO IDFIX
JRST DFIX1 ;GO TO MAIN LINE
HELLO (DFIX,.) ;[235] ENTRY TO DFIX
JRST DFIX1
HELLO (IDINT,.) ;[235] ENTRY TO IDINT
DFIX1:
DMOVE T0,@(L) ;GET THE ARGUMENT
PJRST IDF.0## ;USE IDF.0 ROUTINE
PRGEND
TITLE DBLE %4.(235) REAL TO DOUBLE PRECISION CONVERSION
SUBTTL H. P. WEISS/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1973,1977 BY DIGITAL EQUIPMENT CORPORATION
SEARCH FORPRM
ENTRY DBLE
EXTERN DBLE.
DBLE=DBLE.
PRGEND
TITLE DBLE. %4.(235) REAL TO DOUBLE PRECISION CONVERSION
SUBTTL ED YOURDON/KK
SUBTTL D. TODD/DRT/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1972,1977 BY DIGITAL EQUIPMENT CORPORATION
;FROM 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
;THE CALLING SEQUENCE FOR THE ROUTINE IS AS FOLLOWS:
; JSA Q, DBLE
; EXP ARG
;THE HIGH ORDER PORTION OF THE DOUBLE PRECISION ANSWER IS
;LEFT IN ACCUMULATOR A, AND THE LOW ORDER PORTION IN B.
SEARCH FORPRM
A= 0
B= 1
Q= 16
HELLO (DBLE,.) ;[235] ENTRY TO DOUBLE ROUTINE
MOVE A, @(Q) ;GET THE ARGUMENT IN HIGH ORDER
MOVEI B, 0 ;SET LOW ORDER PORTION TO ZERO
GOODBY (1) ;EXIT
PRGEND
TITLE SNGL %4.(235) DOUBLE PRECISION TO REAL CONVERSION
SUBTTL H. P. WEISS/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1973,1977 BY DIGITAL EQUIPMENT CORPORATION
SEARCH FORPRM
ENTRY SNGL
EXTERN SNGL.
SNGL=SNGL.
PRGEND
TITLE SNGL. %4.(235) DOUBLE PRECISION TO REAL CONVERSION
SUBTTL D. TODD/DMN/TWE/DRT/HPW 11-DEC-73
;THIS SOFTWARE IS FURNISHED UNDER A LICENSE AND MAY ONLY BE USED
; OR COPIED IN ACCORDANCE WITH THE TERMS OF SUCH LICENSE.
;COPYRIGHT (C) 1972,1977 BY DIGITAL EQUIPMENT CORPORATION
;DOUBLE PRECISION TO SINGLE PRECISION CONVERSION FUNCTION
;THIS ROUTINE OBTAINS THE MOST SIGNIFICANT PART OF A DOUBLE
;PRECISION ARGUMENT. THE ROUTINE IS CALLED IN THE FOLLOWING MANNER:
; JSA Q,SNGL
; EXP ARG
;WHERE ARG IS THE ADDRESS OF THE HIGH ORDER PORTION OF A DOUBLE
;PRECISION ARGUMENT, THE LOW ORDER PART BEING IN ARG+1. THE
;ANSWER IS RETURNED IN ACCUMULATOR A.
Q=16
SEARCH FORPRM
HELLO (SNGL,.) ;[235] ENTRY TO SNGL ROUTINE
DMOVE T0,@(L) ;GET THE ARGUMENT IN AC 0
PJRST SNG.0## ;USE AC0 SINGLE ROUTINE
END