; FUNCTION.ASM
; ------------
;
; See FALCONER.WS4 as doc.
;
; (Retyped by Emmanuel ROCHE.)
;
;--------------------------------
; External routines required, see FLTARITH
;--------------------------------
;
; External arithmetic error trap
;
	extrn	aerc
;
; External floating arithmetic
;
	extrn	fadd,fdiv,fdivr
	extrn	fint,fixr
	extrn	fmul,frcip,fsubr
;
; External format conversion
;
	extrn	flot,flota,flotd
;
; External tests and manipulation
;
	extrn	fxchg,fcmp
;
; External memory access
;
	extrn	fload,lfbs
;
;--------------------------------
; Entry points allowed
;
; Functions
;
	entry	fract,fmod,poly
	entry	log2,logb,exp2,expx
;
; Logical operators
;
	entry	.or., .and., .xor.
;
; Relational operators
;
	entry	.equ., .ne., .le., .gt.
	entry	.lt., .ge., .gg.
;
;--------------------------------
; Macro definitions
;--------------------------------
;
; Load (and POP) real (reg) from TOS,
; stored by SFTS macro
;
lfts	macro	reg
bc.l	equ	b
de.h	equ	d
	pop	reg
	mov	5-reg/2,reg
	pop	reg
	if	reg*(reg-d)
	error	"R"
	endif
	endm
;
; Move operation on register pair B, D, or H
;
movd	macro	r1,r2
	if	((r1-d)*(r1-h)*r1) OR ((r2-d)*(r2-h)*r2)
	error	"R"
	endif
	mov	r1,r2
	mov	r1+1,r2+1
	endm
;
; Move floating operand from reg1 to reg2
;
movf	macro	r2,r1
bc.l	equ	b
de.h	equ	d
	if	(r1*(r1-d)) or (r2*(r2-d))
	error	"R"
	endif
	mov	r2,r1
	mov	r2+1,r1+1
	mov	5-r2/2,5-r1/2
	endm
;
; Reload (BC.L), stored by PUSH B, PUSH H sequence
;
reload	macro	reg
bc.l	equ	b
	if	reg-b
	error	'R'
	db	0,0,0
	endif
	if	reg-b=0		; Was IFZ
	pop	b
	mov	l,c
	pop	b
	endif
	endm
;
; "Return" and check stack level zero
;
rtn	macro
	if	.lvl
	error	"0"+.lvl
.lvl	set	0
	endif
	ret
	endm
;
; Save (BC.L), to be restored by RELOAD BC.L later
;
save	macro	reg
bc.l	equ	b
	if	reg-b
	error	"R"
	db	0,0
	endif
	if	reg-b=0		; Was IFZ
	push	b
	push	h
	endif
	endm
;
; Store real value on top of stack;
;   note SFTS B affects (A)
;
sfts	macro	r
bc.l	equ	b
de.h	equ	d
	if	r*(r-d)
	error	"R"
	endif
	if	r=0		; Was IFZ
	mov	a,l
	endif
	push	r
	push	psw-r
	endm
;
;--------------------------------
; Start the code
;--------------------------------
;
; Extract fractional part of (DE.H)
; A,F,D,E,H
;
fract:	save	bc.l		; 
	sfts	d		; 
	call	fint		; 
	lfts	b		; Orig value to (BC.L)
	call	fsubr		; Remove integer portion
	reload	bc.l		; 
	rtn			; 
;
; Convert (DE.H) and (BC.L) to rounded integers
;   in (BC) and (DE) respectively
; A,F,B,C,D,E,H,L
;
fixrt:	call	fixr		; 
	rc			; Overflow
	call	fxchg		; 
	jmp	fixr		; 
;
; Modulo arithmetic
; (DE.H) := (BC.L) modulo (DE.H)
; System trap for (DE.H) = 0
; A,F,D,E,H
;
fmod:	save	bc.l		; 
	sfts	d		; 
	mov	a,b		; 
	xra	d		; Compare signs
	push	psw		; Save for exit
	call	fdivr		; BC.L / DE.H
@01	set	.lvl		; 
	jc	fmod2		; Overflow
	call	fint		; Integer (BC.L / DE.H)
	jc	fmod2		; Overflow
	pop	psw		; 
	jp	fmod1		; Signs same
	lxi	b,8000H		; 
	mvi	l,81H		; -1.000
	call	fadd		; Correct
fmod1:	lfts	b		; Original DE.H
	call	fmul		; DE.H * integer (BC.L / DE.H)
	reload	bc.l		; 
	jnc	fsubr		; BC.L - DE.H * integer (BC.L / DE.H)
	rtn			; FMUL overflowed
.lvl	set	@01		; 
fmod2:	pop	psw		; 
	lfts	d		; Overflow occurred,
	stc			;   restore input condition.
	reload	bc.l		; 
	rtn			; 
;
; Convert (DE.H) to logarithm, base 2
; Trap if (DE.H) <= 0 i.e., error
; Time approx 9 millisec
; A,F,D,E,H
;
log2:	mov	a,h		; 
	ora	a		; 
	cz	aerc		; Zero, trap
	rc			; 
	mov	a,d		; 
	ral			; 
	cc	aerc		; Negative, trap
	rc			; 
	save	bc.l		; 
	push	h		; Save exponent
	movd	b,d		; X to BC.L
	lxi	d,3502H		; SQRT(2)
	lxi	h,8181H		; X range 1 to 2
	call	fadd		; X + SQRT(2)
	sfts	d		; and save
	lxi	d,0B502H	; -SQRT(2)
	mvi	h,81H		; 
	call	fadd		; X-SQRT(2)
	lfts	b		; 
	call	fdiv		; Form term
	movf	b,d		; and copy
	call	fmul		; 
	call	fmul		; Term^3
	sfts	d		; 
	lxi	d,38A6H		; 2.8052
	mvi	h,82H		; 
	call	fmul		; 
	lfts	b		; 
	sfts	d		; 
	lxi	d,7E08H		; 0.9935
	mvi	h,80H		; 
	call	fmul		; 
	lfts	b		; 
	call	fadd		; 
	movf	b,d		; Partial term of BC.L
	pop	d		; Get exponent
	mov	a,d		; 
	sui	81H		; 
	call	flota		; Convert
	call	fadd		; Add characteristic in
	lxi	b,0		; 
	mvi	l,80H		; 0.5000
	call	fadd		; 
	reload	bc.l		; 
	rtn			; 
;
; Log (DE.H) base (BC.L) => (DE.H)
; Carry for overflow. Returns max values, or 0.
; (BC.L) or (DE.H) <= 0 causes trap
; Time approx 20 millisec
; A,F,D,E,H
;
logb:	save	bc.l		; 
	sfts	d		; 
	movf	d,b		; 
	call	log2		; 
	lfts	b		; Restore operand
	sfts	d		; Save log of base
	movf	d,b		; 
	call	log2		; 
	lfts	b		; Restore log base
	call	fdiv		; 
	reload	bc.l		; 
	rtn			; 
;
; Evaluate polynomial in (DE.H) = x
; (DE.H) := A(N)*X^N + A(N-1)*X^(N-1) + ... + A(1)*X + A(0)
; Carry for arithmetic overflow
; (BC) specifies address of coefficients
; First coefficient is order of polynomial (128 max)
; A,F,D,E,H
;
poly:	save	bc.l		; 
	ldax	b		; Get order
	inx	b		; Advance coeff pointer
	sfts	d		; Save argument
@arg	set	.lvl		; Argument stack address
	mvi	h,0		; Clear partial value
	push	psw		; Save order counter
poly1:	push	b		; Save coeff loc
	sfts	d		; Save partial value
	call	fload		; Get coefficient
	lfts	b		; Recover partial value to (BC.H)
	call	fadd		; Add in
	pop	b		; Coeff pointer
	jc	poly2		; Arith overflow
	pop	psw		; Order counter
	dcr	a		; 
	jm	poly3		; Done
	push	psw		; Save order counter
	push	b		; Save coeff pointer
	mvi	a,.lvl-@arg	; 
	call	lfbs		; Get argument
	call	fmul		; Multiply
	pop	b		; Restore coeff pointer
	inx	b		; 
	inx	b		; 
	inx	b		; Advance to next coeff
	jnc	poly1		; No arith error
poly2:	pop	b		; Error exit, purge stack
poly3:	pop	b		; 
	pop	b		; Purge argument from stack
	reload	bc.l		; 
	rtn			; 
;
; Exponential (DE.H) := 2^(DE.H)
; Carry for overflow
; A,F,D,E,H
;
exp2:	mov	a,d		; 
	ora	a		; 
	jp	exp21		; 
	xri	80H		; Set positive
	mov	d,a		; 
	call	exp21		; 
	cnc	frcip		; Neg exponent
	rnc			; 
	mvi	h,0		; Zero for negative overflow
	rtn			; 
exp21:	save	bc.l		; 
	movf	b,d		; Copy argument to B
	call	fixr		; 
	jc	exp22		; Too large, overflow
	push	d		; Save integer portion
	call	flotd		; 
	call	fsubr		; Form fractional portion
	lxi	b,ex2c		; Point to coefficients
	call	poly		; Form 2^(fract(x))
	movf	b,d		; 
	call	fmul		; Form (1+A1*X+...+AN*X^N)^2
	pop	b		; Get integer portion(x)
	mov	a,b		; 
	ora	a		; 
	stc			; 
	jnz	exp22		; Too large, overflow
	mov	a,c		; 
	add	h		; 
	mov	h,a		; Exponent overlow causes Carry
exp22:	reload	bc.l		; 
	rnc			; 
	lxi	d,7FFFH		; 
	mov	h,e		; Set max value
	rtn			; 
;
; Polynomial coefficients for 2^(x)
;
ex2c:	db	3		; Polynomial order
	db	7AH,01H,06H	; 0.0081790
	db	7CH,0DH,73H	; 0.059340
	db	7FH,81H,31H	; 0.34669
	db	81H,00H,00H	; 1.0000
;
; Exponential (DE.H) := (BC.L)^(DE.H)
; (BC.L) < 0 illegal, divertto trap.
; (BC.L) and (DE.H) = 0 illegal, trap.
; Carry for over/underflow, returns max, 0.1
; A,F,D,E,H
;
expx:	mov	a,l		; 
	ora	a		; 
	jnz	expx1		; (BC.L) <> 0
	ora	h		; 
	cz	aerc		; Illegal, trap
	mvi	h,0		; 0^any = 0
	rtn			; 
expx1:	mov	a,b		; 
	ora	a		; 
	cm	aerc		; Illegal, trap
	rc			; 
	mov	a,h		; 
	ora	a		; 
	jnz	expx3		; 
expx2:	lxi	d,0		; 
	mvi	h,81H		; Any^0 = 1.000
	rtn			; 
expx3:	save	bc.l		; 
	sfts	d		; 
	movf	d,b		; 
	call	log2		; 
	lfts	b		; Restore argument
	call	fmul		; 
	reload	bc.l		; 
	jnc	exp2		; 
	mov	a,h		; 
	ora	a		; 
	stc			; 
	jz	expx2		; Underflow, return 1.000
	mov	a,d		; 
	ora	a		; 
	stc			; 
	rp			; +ve overflow, return max
	mvi	h,0		; -ve overflow, return 0
	rtn			; 
;
;--------------------------------
; The logical operators
;   treat all arguments as signed integers
;   and return the floating representation of
;   the bitwise operation specified.
;   Error: If any argument is outside the
;   range -32768 to 32767.
;--------------------------------
;
; Logical OR on (BC.L),(DE.H)
; A,F,D,E,H
;
.or.:	save	bc.l		; 
	call	fixrt		; 
	jc	.or.2		; 
	mov	a,b		; 
	ora	d		; 
	mov	d,a		; 
	mov	a,c		; 
	ora	e		; 
.or.1:	mov	e,a		; 
	call	flot		; 
.or.2:	reload	bc.l		; 
	rtn			; 
;
; Logical AND on (BC.L),(DE.H)
; A,F,D,E,H
;
.and.:	save	bc.l		; 
	call	fixrt		; 
	jc	.or.2		; 
	mov	a,b		; 
	ana	d		; 
	mov	d,a		; 
	mov	a,c		; 
	ana	e		; 
	jmp	.or.1		; 
.lvl	set	.lvl-2		; 
;
; Logical XOR on (BC.L),(DE.H)
; A,F,D,E,H
;
.xor.:	save	bc.l		;
	call	fixrt		; 
	jc	.or.2		; 
	mov	a,b		; 
	xra	d		; 
	mov	d,a		; 
	mov	a,c		; 
	xra	e		; 
	jmp	.or.1		; 
.lvl	set	.lvl-2		; 
;
;--------------------------------
; The relational operators
;   return -1 for true
;           0 for false.
;--------------------------------
;
; Test (DE.H) = (BC.L)
; If so, (DE.H) := -1, else 0
; A,F,D,E,H
;
.equ.:	call	fcmp		; 
	jz	.tru		; True
.mtru:	xra	a		; 
	mov	h,a		; (DE.H) := 0
	rtn			; 
;
; Test (DE.H) <> (BC.L)
; If so, (DE.H) := -1, else 0
; A,F,D,E,H
;
.ne.:	call	fcmp		; 
	jz	.mtru		; False
.tru:	mvi	h,81H		; (DE.H) := -1.0
	lxi	d,8000H		; Use LXI D,0 for true = +1.0, for Pascal etc
	ora	a		; Clear any Carry
	rtn			; 
;
; Test (DE.H) <= (BC.L)
; If so, (DE.H) := -1, else 0
; A,F,D,E,H
;
.le.:	call	fcmp		; 
.le.1:	jp	.tru		; 
	jmp	.mtru		; 
;
; Test (DE.H) > (BC.L)
; If so, (DE.H) := -1, else 0
; A,F,D,E,H
;
.gt.:	call	fcmp		; 
.gt.1:	jm	.tru		; 
	jmp	.mtru		; 
;
; Test (DE.H) < (BC.L)
; If so, (DE.H) := -1, else 0
; A,F,D,E,H
;
.lt.:	call	fcmp		; 
	jz	.mtru		; 
	jmp	.le.1		; 
;
; Test (DE.H) >= (BC.L)
; If so, (DE.H) := -1, else 0
; A,F,D,E,H
;
.ge.:	call	fcmp		; 
	jz	.tru		; 
	jmp	.gt.1		; 
;
; Set (DE.H) := (DE.H) * 2^15 and perform .GT.
;
; This can be used to test for a value effectively
; zero with respect to another. For termination of
; iteration loops, etc.  The value 2^15 applies to
; this arithmetic system, and should be customized
; to the precision of any particular arithmetic
; system for program portability.
;
.gg.:	mov	a,h		; 
	adi	15		; 
	jnc	.gg.1		; Dynamic room
	push	h		; 
	mov	a,l		; 
	sui	15		; 
	mov	l,a		; 
	jnc	.gg.2		; Dynamic room
	pop	h		; No room, return false
	jmp	.mtru		; 
.gg.1:	mov	h,a		; 
	push	h		; 
.gg.2:	call	fcmp		; 
	pop	h		; 
	jmp	.gt.1		; 
;
;--------------------------------
;
	end			; of FUNCTION.ASM