home *** CD-ROM | disk | FTP | other *** search
- *
- * Sqrt.asm (of PCQ Pascal)
- * Copyright 1990 Patrick Quaid
- *
- * This routine just implements the sqrt function using
- * Newton's method. As I recall, this routine will in
- * some cases diverge, but it seems to work reasonably
- * well. It figures the square root to 5 or 6 digits,
- * which in this implementation requires between, say,
- * 6 and 12 iterations (for 100 and 200000, respectively).
- *
-
- XREF _p%MathBase
-
- XREF _LVOSPDiv
- XREF _LVOSPMul
- XREF _LVOSPCmp
- XREF _LVOSPAdd
- XREF _LVOSPSub
- XREF _LVOSPAbs
-
- XDEF _p%sqrt
- _p%sqrt
- move.l 4(sp),d0 ; d0 := x
- move.l #$80000042,d1 ; d1 := 2.0
- move.l _p%MathBase,a6
- jsr _LVOSPDiv(a6) ; d0 := x / 2.0
- move.l d0,-(sp) ; save approx on the stack
-
- move.l 8(sp),d0 ; d0 := x
- move.l #$C3500051,d1 ; d1 := 100000.0
- jsr _LVOSPDiv(a6) ; d0 := x / 100000.0
- move.l d0,-(sp) ; save epsilon on the stack
-
- ComputeError:
- move.l 4(sp),d0 ; get approx
- move.l d0,d1 ; move to d1
- jsr _LVOSPMul(a6) ; d0 := Sqr(approx)
- move.l d0,d1 ; d1 := Sqr(approx)
- move.l 12(sp),d1 ; d0 := x
- jsr _LVOSPSub(a6) ; d0 := x - Sqr...
- jsr _LVOSPAbs(a6) ; d0 := abs(x - sqr())
- move.l (sp),d1 ; d1 := epsilon
- jsr _LVOSPCmp(a6)
- ble.s ReturnApprox ; if error <= epsilon, leave
-
- move.l 12(sp),d0 ; d0 := x
- move.l 4(sp),d1 ; d1 := approx
- jsr _LVOSPDiv(a6) ; d0 := x / approx
- move.l 4(sp),d1 ; d1 := approx
- jsr _LVOSPAdd(a6) ; d0 := approx + x / approx
- move.l #$80000042,d1 ; d1 := 2.0
- jsr _LVOSPDiv(a6) ; d0 := (approx + x/approx)/2
- move.l d0,4(sp) ; approx := d0
- bra ComputeError ; iterate
- ReturnApprox
- move.l 4(sp),d0 ; d0 := approx
- addq.l #8,sp ; pop stack
- rts
-
- END
-