Web pdp-10.trailing-edge.com

Trailing-Edge - PDP-10 Archives - decuslib10-10 - 43,50517/sqrt.mac
There is 1 other file named sqrt.mac in the archive. Click here to see a list.
```	TITLE	SQRT for RPGLIB V1
SUBTTL	Square root and floating point routines

;	SQRT for RPGLIB V1
;
;
;	This module contains the square root routine along with other
;	miscellaneous floating point routines.
;
;	Copyright (C) 1976, Bob Currier and Cerritos College
;
;

SEARCH	RPGPRM, UUOSYM, MACTEN

%%LBLP==:%%LBLP
DEBUG==:DEBUG
BIS==:BIS

SALL

TWOSEG
RELOC	400000

ENTRY	SQRT.			; square root routine
ENTRY	FLOT.1			; float a single precision integer
ENTRY	FLOT.2			; float a double precision integer
ENTRY	FIX.			; fix a double precision integer

PR==15
NM==5
AC0==0
AC1==1
G==T2
H==T3
```
```;SQRT.		Main square root routine
;
;Uses Newton's method
;
;

SQRT.:	JUMPL	AC1,SQRT.1		; must be positive
MOVE	G,AC1			; get arg
FDVR	G,[2.0]			; get a guess

SQ.1:	MOVE	H,AC1			; get arg
FDVR	H,G			; get H/G
FDVR	H,[2.0]			; get (G+(H/G))/2
MOVE	T1,H			; get it where we can play with it
FSB	T1,G			; get difference
TLNN	T1,377000		; get it to nearest 1e-127
JRST	SQ.2			; yes -
MOVE	G,H			; no get new guess
JRST	SQ.1			; and loop

SQ.2:	MOVE	AC1,H			; get result into proper AC
POPJ	PP,			; and exit

SQRT.1:	PUSHJ	PP,%%H.11##		; square root of negative number illegal
SETZ	AC1,			; return 0 on continue
POPJ	PP,			; exit
```
```;FIX.		Convert floating point to double precision fixed point
;
;Calling sequence:
;
;	MOVE	AC,[floating.point.number]
;	PUSHJ	17,FIX.
;
;

FIX.:	LDB	T1,[POINT 4,PARM,12]	; get the AC
MOVE	T1,(T1)			; get the floating point number
TDNN	T1,[XWD 777,777777]	; is it zero?
JRST	CC3C2D			; yes -
MOVM	T2,T1			; no -
MULI	T2,400			; get exponent
MOVEI	T4,0			; zap
ASHC	T3,-306(T2)		; shift to get integer
JUMPGE	T1,CC3C2C		; jump if input was positive
SETCA	T3,			; else negate the two word result
MOVNS	T4			; complement and increment
SKIPN	T4			; zero?
ADDI	T3,1			; yes - increment
SKIPGE	T3
TLO	T4,1B18			; set sign bit if necessary

CC3C2C:	MOVEM	T3,0(PARM)		; stash result
MOVEM	T4,1(PARM)		; both parts
POPJ	PP,			; exit

CC3C2D:	SETZM	0(PARM)			; input is zero - give zero result
SETZM	1(PARM)			; double precision
POPJ	PP,			; exit
```
```;FLOT.1		Convert single precision fixed point to floating point
;
;Calling sequence:
;
;	MOVE	16,[Z AC,operand]
;	PUSHJ	17,FLOT.1
;
;Exit with result in accumulator specified by AC field of 16
;
;

FLOT.1:	SKIPN	NM,0(PARM)		; is input zero?
JRST	FLOT1A			; yes -
IDIVI	NM,1B18			; no - break it into two pieces
CAIE	NM,0			; high part zero?
TLC	NM,254000		; no - jam exponent
TLC	NM+1,233000		; jam exponent of low part

FLOT1A:	LDB	T2,[POINT 4,PARM,12]	; get AC field
MOVEM	NM,0(T2)		; store result
POPJ	PP,			; and exit
```
```;FLOT.2		Convert double precision fixed point number to floating point
;
;Calling sequence:
;
;	PUSHJ	17,FLOT.2
;
;Exit with result in AC specified by AC field of 16
;
;

FLOT.2:	MOVE	PR,[1.0]		; set scaling factor
MOVE	T2,1(PARM)		; get low order
SKIPL	T1,0(PARM)		; get high order part - is it negative?
JRST	FLOT2B			; no -
TLO	T2,1B18			; yes - be sure low part is negative also
SETCA	T1,			; negate the two word value
MOVMS	T2			; get magnitude
TLZ	T2,1B18			; reset sign bit
JUMPN	T2,FLOT2C		; zero?
AOJA	T1,FLOT2C		; yes - bump and jump

FLOT2B:	TLZ	T2,1B18			; make sure low part is positive

FLOT2C:	TLNE	T1,377776		; is value > 2**54 ?
JRST	FLOT2E			; yes -
ASHC	T1,10			; no - shift over 8 bits
ASH	T2,-10			; restore low part
SKIPE	T1			; set exponent in high part unless
TLC	T1,266000		;   it is zero
SKIPE	T2			; set exponent in low part unless
TLC	T2,233000		;   it is negative
FMP	T1,PR			; multiply by scaling factor
SKIPGE	0(PARM)			; was input negative?
MOVNS	T1			; yes - negate result
LDB	T2,[POINT 4,PARM,12]	; get AC field
MOVEM	T1,0(T2)		; store result
POPJ	PP,			; exit

FLOT2E:	ADD	PR,[1B8]		; number > 2**54 - increase scaling factor
ASHC	T1,-1			; reduce T1
JRST	FLOT2C			; try again
```
```

END
```