ITS Muddle.
[pdp10-muddle.git] / MUDDLE / arith.58
diff --git a/MUDDLE/arith.58 b/MUDDLE/arith.58
new file mode 100644 (file)
index 0000000..1e1e933
--- /dev/null
@@ -0,0 +1,626 @@
+TITLE ARITHMETIC PRIMITIVES FOR MUDDLE
+
+;BKD
+
+;DEFINES MUDDLE PRIMITIVES:   FIX,FLOAT,ATAN,IEXP,LOG,
+;      G?,L?,0?,1?,+,-,*,/,MAX,MIN,ABS,SIN,COS,SQRT,RANDOM,
+;      TIME,SORT.
+
+RELOCATABLE
+
+.INSRT MUDDLE >
+
+O=0
+
+
+DEFINE TYP1
+       (AB) TERMIN
+DEFINE VAL1
+       (AB)+1 TERMIN
+
+DEFINE TYP2
+       (AB)+2 TERMIN
+DEFINE VAL2
+       (AB)+3 TERMIN
+
+DEFINE TYP3
+       (AB)+4 TERMIN
+DEFINE VAL3
+       (AB)+5 TERMIN
+
+DEFINE TYPN
+       (D) TERMIN
+DEFINE VALN
+       (D)+1 TERMIN
+
+
+YES:   MOVSI   A,TATOM ;RETURN PATH FOR 'TRUE'
+       MOVE    B,MQUOTE T
+       JRST FINIS
+
+NO:    MOVSI   A,TFALSE        ;RETURN PATH FOR 'FALSE'
+       MOVEI   B,NIL
+       JRST FINIS
+
+\f;ERROR RETURNS AND OTHER UTILITY ROUTINES
+
+OVRFLW==10
+OVRFLD:        PUSH    TP,$TATOM
+       PUSH    TP,MQUOTE OVERFLOW
+       JRST    CALER1
+
+ARGCHK:                        ;CHECK FOR SINGLE FIXED OR FLOATING
+                       ;ARGUMENT IF FIXED CONVERT TO FLOATING
+                       ;RETURN FLOATING ARGRUMENT IN B ALWAYS
+       ENTRY   1
+       HLRZ    C,TYP1  
+       MOVE    B,VAL1
+       CAIN    C,TFLOAT        ;FLOATING?
+       POPJ P, ;YES, RETURN
+       CAIE    C,TFIX  ;FIXED?
+       JRST    WTYP    ;NO, ERROR
+       JSP A,BFLOAT    ;YES, CONVERT TO FLOATING AND RETURN
+       POPJ P,
+
+OUTRNG:        PUSH    TP,$TATOM
+       PUSH    TP,MQUOTE ARGUMENT-OUT-OF-RANGE
+       JRST    CALER1
+
+NSQRT: PUSH TP,$TATOM
+       PUSH TP,MQUOTE NEGATIVE-ARGUMENT
+       JRST CALER1"
+
+WTYP:  PUSH TP,$TATOM
+       PUSH TP,MQUOTE WRONG-TYPE
+       JRST CALER1
+
+DEFINE MFLOAT AC
+       IDIVI AC, 400000
+       FSC     AC+1,233
+       FSC     AC,254
+       FADR AC,AC+1
+       TERMIN
+
+BFLOAT:        MFLOAT  B
+       JRST    (A)
+
+OFLOAT:        MFLOAT  O
+       JRST    (C)
+
+BFIX:  MULI    B,400
+       TSC     B,B
+       ASH     C,(B)-243
+       MOVE    B,C
+       JRST    (A)
+
+\f;DISPATCH TABLES USED TO CONTROL THE FLOW OF THE VARIOUS PRIMITIVES
+
+TABLE2:        NO      ;TABLE2 (0)
+TABLE3:        YES     ;TABLE2 (1)  &  TABLE3 (0)
+       NO      ;TABLE2 (2)
+
+
+FUNC:          JSP     A,BFIX
+       JSP     A,BFLOAT
+       SUB     B,VALN
+       IDIV    B,VALN
+       ADD     B,VALN
+       IMUL    B,VALN
+       JSP     C,SWITCH
+       JSP     C,SWITCH
+
+FLFUNC==.-2
+       FSBR    B,O
+       FDVR    B,O
+       FADR    B,O
+       FMPR    B,O
+       JSP     C,FLSWCH
+       JSP     C,FLSWCH
+\f;PRIMITIVES FLOAT AND FIX
+
+MFUNCTION      FIX,SUBR
+       MOVEI   E,0
+       JRST    TRANS
+
+MFUNCTION      FLOAT,SUBR
+       MOVEI   E,1
+
+TRANS: ENTRY   1
+       MOVE    A,TYP1
+       MOVE    B,VAL1
+       CAMN    A,TYPS(E)+1     ;SAME TYPE ARGUMENT?
+       JRST    FINIS
+       CAME    A,TYPS(E)       ;correct type argument ?
+       JRST    WTYP
+       XCT     FUNC(E) ;perform appropriate operation
+       MOVE    A,TYPS(E)+1     ;save this new type
+JRST FINIS
+
+TYPS:          TFLOAT,,0
+       TFIX,,0
+       TFLOAT,,0
+
+MFUNCTION      ABS,SUBR
+       ENTRY   1
+       MOVE    A,TYP1
+       CAME    A,$TFIX
+       CAMN    A,$TFLOAT
+       JRST    MOVIT
+       JRST    WTYP
+MOVIT: MOVM    B,VAL1  ;GET ABSOLUTE VALUE OF ARGUMENT
+       JRST    FINIS
+
+MFUNCTION      MOD,SUBR
+       ENTRY   2
+       MOVSI   A,TFIX
+       CAME    A,TYP1  ;FIRST ARG FIXED ?
+       JRST    WTYP
+       CAME    A,TYP2  ;SECOND ARG FIXED ?
+       JRST    WTYP
+       MOVE    B,VAL1
+       IDIV    B,VAL2  ;FORM QUOTIENT & REMAINDER
+       JUMPGE  C,.+2   ;Only return positive remainders
+       ADD     C,VAL2
+       MOVE    B,C     ;RETURN REMAINDER
+       JRST    FINIS
+\f;PRIMITIVES PLUS, DIFFERENCE, TIMES, DIVIDE, MIN, AND MAX
+
+MFUNCTION      MIN,SUBR
+       MOVEI   E,6
+       JRST    GOPT
+
+       MFUNCTION       MAX,SUBR
+       MOVEI   E,7
+GOPT:  ENTRY
+       MOVE    D,AB    ;ARGUMENT POINTER
+       JUMPL   D,MINMAX        ;ANY ARGUMENTS AT ALL ?
+       MOVSI   A,TFLOAT        ;DEFAULT TYPE
+       MOVE    B,INFIN(E)      ;DEFAULT VALUE + OR - "LARGE NUMBER"
+       JRST    FINIS
+INFIN==.-6
+       377777,,-1
+       400000,,1
+
+MFUNCTION      DIVIDE,SUBR,[/]
+       MOVEI   E,3
+       JRST    ARITH0
+
+MFUNCTION      DIFFERENCE,SUBR,[-]
+       MOVEI   E,2
+       JRST    ARITH0
+
+MFUNCTION      TIMES,SUBR,[*]
+       MOVEI   E,5
+       JRST    ARITH0
+
+MFUNCTION      PLUS,SUBR,[+]
+       MOVEI   E,4
+
+ARITH0:        ENTRY
+       MOVE    D,AB    ;argument pointer
+       CAMGE   D,[-2,,0]       ;LESS THAN TWO ARGUMENTS ?
+       JRST    MINMAX
+       MOVSI   A,TFIX  ;initial type of result
+       MOVE    B,E     ;initial accumulator contents for zero & one argument
+       TRZ     B,-2
+       JRST    MINMAX+3
+MINMAX:        MOVE    A,TYP1
+       MOVE    B,VAL1  ;initial value of accumulator for more than one argument is first value
+       ADD     D,[2,,2]        ;UPDATE ARGUMENT POINTER
+       JUMPGE  D,FINIS ;ANY MORE ARGUMENTS ?
+       JFCL    OVRFLW,.+1
+       CAME    A,$TFIX ;WAS THE FIRST ARGUMENT FIXED ?
+       JRST    ARITH3
+ARITH1:        CAME    A,TYPN  ;next argument fixed ?
+       JRST    ARITH2
+       XCT     FUNC(E) ;PERFORM APPROPRIATE OPERATION
+       ADD     D,[2,,2]        ;UPDATE ARGUMENT POINTER
+       JUMPL   D,ARITH1        ;repeat for next argument if any
+       JFCL    OVRFLW,OVRFLD
+       JRST    FINIS
+\f;CONTINUATION OF PLUS,TIMES, ETC.
+
+ARITH3:        CAME    A,$TFLOAT       ;was the first argument floating ?
+       JRST    WTYP
+       SKIPA
+
+ARITH2:        JSP     A,BFLOAT        ;float accumulator contents
+       MOVE    C,TYPN  ;get next argument's type
+       MOVE    O,VALN  ;get next argument's value
+       CAMN    C,$TFLOAT       ;floating ?
+       JRST    OPERATE
+       CAME    C,$TFIX ;fixed ?
+       JRST    WTYP
+       JSP     C,OFLOAT        ;go float this fixed argument
+OPERATE:       XCT     FLFUNC(E)       ;perform appropriate operation
+       ADD     D,[2,,2]        ;UPDATE ARGUMENT POINTER
+       JUMPL   D,ARITH2+1      ;repeat for next argument if any
+       JFCL    OVRFLW,OVRFLD
+       MOVSI   A,TFLOAT
+       JRST FINIS
+
+SWITCH:        XCT     COMPAR(E)       ;FOR MAX & MIN TESTING
+       MOVE    B,VALN
+       JRST    (C)
+COMPAR==.-6
+       CAMLE   B,VALN
+       CAMGE   B,VALN
+
+FLSWCH:        XCT     FLCMPR(E)
+       MOVE    B,O
+       JRST    (C)
+FLCMPR==.-6
+       CAMLE   B,O
+       CAMGE   B,O
+\f;PRIMITIVES ONEP AND ZEROP
+
+MFUNCTION      ONEP,SUBR,[1?]
+       MOVEI   E,1
+       JRST    JOIN
+
+MFUNCTION      ZEROP,SUBR,[0?]
+       MOVEI   E,
+
+JOIN:  ENTRY 1
+       MOVE    A,TYP1
+       CAMN    A,$TFIX ;fixed ?
+       JRST    TESTFX
+       CAME    A,$TFLOAT       ;floating ?
+       JRST    WTYP
+       MOVE    B,VAL1
+       CAMN    B,NUMBR(E)      ;equal to correct value ?
+       JRST    YES
+       JRST    NO
+
+TESTFX:        CAMN    E,VAL1  ;equal to correct value ?
+       JRST    YES
+       JRST    NO
+
+NUMBR: 0       ;FLOATING PT  ZERO
+       201400,,0       ;FLOATING PT ONE
+\f;PRIMITIVES LESSP AND GREATERP
+
+
+MFUNCTION      LESSP,SUBR,[L?]
+       MOVEI   E,1
+       JRST    ARGS
+
+MFUNCTION      GREATERP,SUBR,[G?]
+       MOVEI   E,0
+
+ARGS:  ENTRY 2
+       MOVE    O,VAL1
+       MOVE    A,TYP1
+       MOVE    B,VAL2
+       SETO    D,      ;used for flow of control in this routine
+       CAMN    A,$TFLOAT
+       AOJA    D,CONT
+       CAME    A,$TFIX
+       JUMPL   D,WTYP
+CONT:  MOVE    A,TYP2
+       CAMN    A,$TFIX
+       AOJE    D,FIXFIX        ;are both arguments fixed
+       CAME    A,$TFLOAT
+       JRST    FLTFIX
+       JUMPE   D,FLTFLT        ;are both arguments floating ?
+       JSP     C,OFLOAT        ;go float the first argument
+FLTFLT:        FSBR    O,B     ;both arguments are floating here
+TEST:  JUMPL   O,@TABLE2(E)
+       JUMPG   O,@TABLE3(E)
+       JRST    NO
+
+FLTFIX:        JUMPLE  D,WTYP
+       JSP     A,BFLOAT        ;go float the second argument
+       JRST    FLTFLT
+
+FIXFIX:        SUB     O,B     ;both arguments are fixed here
+       JRST    TEST
+
+MFUNCTION RANDOM,SUBR
+       ENTRY
+       HLRE    A,AB
+       CAMGE   A,[-4]  ;At most two arguments to random to set seeds
+       JRST    WNA
+       JRST    RANDGO(A)
+       MOVE    B,VAL2  ;Set second seed
+       MOVEM   B,RLOW
+       MOVE    A,VAL1  ;Set first seed
+       MOVEM   A,RHI
+RANDGO:        MOVE B,RLOW     ;FREDKIN'S RANDOM NUMBER GENERATOR.
+       MOVE A,RHI
+       MOVEM A,RLOW
+       LSHC A,-43
+       XORB B,RHI
+       MOVSI A,TFIX
+       JRST FINIS
+RHI:   267762113337
+RLOW:  155256071112
+\fMFUNCTION SQRT,SUBR
+       ENTRY 1
+       MOVE B,1(AB)
+       HLRZ A,(AB)
+       CAIN A,TFLOAT
+       JRST SQ1
+       CAIE A,TFIX
+       JRST WTYP
+       JSP A,BFLOAT
+SQ1:   JUMPL B,NSQRT
+
+       MOVE A,B
+       ASH B,-1
+       FSC B,100
+SQ2:   MOVE C,B        ;NEWTON'S METHOD, SPECINER'S HACK.
+       FDVRM A,B
+       FADRM C,B
+       FSC B,-1
+       CAME C,B
+       JRST SQ2
+       MOVSI A,TFLOAT
+       JRST FINIS
+
+
+MFUNCTION COS,SUBR
+       ENTRY 1
+       MOVE B,1(AB)
+       HLRZ A,(AB)
+       CAIN A,TFLOAT
+       JRST COS1
+       CAIE A,TFIX
+       JRST WTYP
+       JSP A,BFLOAT
+COS1:  FADR B,[1.570796326]    ;COS(X)=SIN (X+PI/2)
+       PUSHJ P,.SIN
+       MOVSI A,TFLOAT
+       JRST FINIS
+
+MFUNCTION SIN,SUBR
+       ENTRY 1
+       MOVE B,1(AB)
+       HLRZ A,(AB)
+       CAIN A,TFLOAT
+       JRST SIN1
+       CAIE A,TFIX
+       JRST WTYP
+       JSP A,BFLOAT
+SIN1:  PUSHJ P,.SIN
+       MOVSI A,TFLOAT
+       JRST FINIS
+
+.SIN:  MOVM A,B
+       CAMG A,[.0001]
+       POPJ P,         ;GOSPER'S RECURSIVE SIN.
+       FDVR B,[-3.0]   ;SIN(X)=4*SIN(X/-3)**3-3*SIN(X/-3)
+       PUSHJ P,.SIN
+       FSC A,1
+       FMPR A,A
+       FADR A,[-3.0]
+       FMPRB A,B
+       POPJ P,
+MFUNCTION      LOG,SUBR
+       PUSHJ P,ARGCHK  ;LEAVES ARGUMENT IN B
+       JUMPLE  B,OUTRNG
+       LDB     D,[331100,,B]   ;GRAB EXPONENT
+       SUBI    D,201   ;REMOVE BIAS
+       TLZ     B,777000        ;SET EXPONENT
+       TLO     B,201000        ; TO 1
+       MOVE    A,B
+       FSBR    A,RT2
+       FADR    B,RT2
+       FDVB    A,B
+       FMPR    B,B
+       MOVE    C,[0.434259751]
+       FMPR    C,B
+       FADR    C,[0.576584342]
+       FMPR    C,B
+       FADR    C,[0.961800762]
+       FMPR    C,B
+       FADR    C,[2.88539007]
+       FMPR    C,A
+       FADR    C,[0.5]
+
+       MOVE    B,D
+       FSC     B,233
+       FADR    B,C
+       FMPR    B,[0.693147180] ;LOG E OF 2
+       MOVSI   A,TFLOAT
+       JRST    FINIS
+RT2:   1.41421356
+\fMFUNCTION     ATAN,SUBR
+       PUSHJ P,ARGCHK
+       MOVM    D,B
+       CAMG    D,[0.4^-8]      ;SMALL ENOUGH SO ATAN(X)=X?
+       JRST    ATAN3   ;YES
+       CAML    D,[7.0^7]       ;LARGE ENOUGH SO THAT ATAN(X)=PI/2?
+       JRST    ATAN1   ;YES
+       MOVN    C,[1.0]
+       CAMLE   D,[1.0] ;IS ABS(X)<1.0?
+       FDVM    C,D     ;NO,SCALE IT DOWN
+       MOVE    B,D
+       FMPR    B,B
+       MOVE    C,[1.44863154]
+       FADR    C,B
+       MOVE    A,[-0.264768620]
+       FDVM    A,C
+       FADR    C,B
+       FADR    C,[3.31633543]
+       MOVE    A,[-7.10676005]
+       FDVM    A,C
+       FADR    C,B
+       FADR    C,[6.76213924]
+       MOVE    B,[3.70925626]
+       FDVR    B,C
+       FADR    B,[0.174655439]
+       FMPR    B,D     ;
+       JUMPG   D,ATAN2 ;WAS ARG SCALED?
+       FADR    B,PI2   ;YES,  ATAN(X)=PI/2-ATAN(1/X)
+       JRST    ATAN2
+ATAN1: MOVE    B,PI2
+ATAN2: SKIPGE  1(AB)   ;WAS INPUT NEGATIVE?
+       MOVNS   B               ;YES,COMPLEMENT
+ATAN3: MOVSI   A,TFLOAT        
+       JRST    FINIS
+PI2:   1.57079632
+\fMFUNCTION     IEXP,SUBR,[EXP] 
+       PUSHJ P,ARGCHK  ;LEAVE FLOATING POINT ARG IN B
+       MOVM    A,B
+       SETZM   B
+       FMPR    A,[0.434294481] ;LOG BASE 10 OF E
+       MOVE    D,[1.0]
+       CAMG    A,D
+       JRST    RATEX
+       MULI    A,400
+       ASHC    B,-243(A)
+       CAILE   B,43
+       JRST    OUTRNG
+       CAILE   B,7
+       JRST    EXPR2
+EXPR1: FMPR    D,FLOAP1(B)
+       LDB     A,[103300,,C]   
+       SKIPE   A
+       TLO     A,177000
+       FADR    A,A
+RATEX: MOVEI   B,7
+       SETZM   C
+RATEY: FADR    C,COEF2-1(B)
+       FMPR    C,A
+       SOJN    B,RATEY
+       FADR    C,[1.0] 
+       FMPR    C,C
+       FMPR    D,C
+       MOVE    B,[1.0]
+       SKIPL   1(AB)   ;SKIP IF INPUT NEGATIVE
+       SKIPN   B,D
+       FDVR    B,D
+       MOVSI   A,TFLOAT
+       JRST    FINIS
+EXPR2: LDB     E,[030300,,B]   
+       ANDI    B,7
+       MOVE    D,FLOAP1(E)
+       FMPR    D,D     ;TO THE 8TH POWER
+       FMPR    D,D
+       FMPR    D,D
+       JRST    EXPR1
+
+COEF2: 1.15129278
+       0.662730884
+       0.254393575
+       0.0729517367
+       0.0174211199
+       2.55491796^-3
+       9.3264267^-4
+
+FLOAP1:        1.0
+       10.0
+       100.0
+       1000.0
+       10000.0
+       100000.0
+       1000000.0
+       10000000.0
+\f;routine to sort lists or vectors of either fixed point or floating numbers
+;the components are interchanged repeatedly to acheive the sort
+;first arg:    the structure to be sorted
+;if no second arg sort in descending order
+;second arg:   if false then sort in ascending order
+;              else sort in descending order
+
+MFUNCTION      SORT,SUBR
+       ENTRY 
+       HLRZ    A,AB
+       CAIGE   A,-4    ;Only two arguments allowed
+       JRST    WNA
+       MOVE    O,DESCEND       ;Set up "O" to test for descending order as default condition
+       CAIE    A,-4    ;Optional second argument?
+       JRST    .+4
+       HLRZ    B,TYP2  ;See if it is other than false
+       CAIN    B,TFALSE
+       MOVE    O,ASCEND        ;Set up "O" to test for ascending order
+       HLRZ    A,TYP1  ;CHECK TYPE OF FIRST ARGUMENT
+       CAIN    A,TLIST
+       JRST    LSORT
+       CAIN    A,TVEC
+       JRST    VSORT
+       JRST    WTYP
+
+
+
+
+GOBACK:        MOVE    A,TYP1  ;RETURN THE SORTED ARGUMENT AS VALUE
+       MOVE    B,VAL1
+       JRST    FINIS
+
+DESCEND:       CAMG    C,(A)+1
+ASCEND:                CAML    C,(A)+1
+\f;ROUTINE TO SORT LISTS IN NUMERICAL ORDER
+
+LSORT: MOVE    A,VAL1
+       JUMPE   A,GOBACK        ;EMPTY LIST?
+       HLRZ    B,(A)   ;TYPE OF FIRST COMPONENT
+       CAIE    B,TFIX
+       CAIN    B,TFLOAT
+       SKIPA
+       JRST    WTYP
+       MOVEI   E,0     ;FOR COUNT OF LENGTH OF LIST
+LCOUNT:        JUMPE   A,LLSORT        ;REACHED END OF LIST?
+       MOVE    A,(A)   ;NEXT COMPONENT
+       TLZ     A,(B)   ;SAME TYPE AS FIRST COMPONENT?
+       TLNE    A,-1
+       JRST    WTYP
+       AOJA    E,LCOUNT        ;INCREMENT COUNT AND CONTINUE
+
+LLSORT:        SOJE    E,GOBACK        ;FINISHED WITH SORTING?
+       HRRZ    A,VAL1  ;START THIS LOOP OF SORTING AT THE BEGINNING
+       MOVEM   E,(P)+1 ;Save the iteration depth
+CLSORT:        HRRZ    B,(A)   ;NEXT COMPONENT
+       MOVE    C,(B)+1 ;ITS VALUE
+       XCT     O       ;ARE THESE TWO COMPONENTS IN ORDER?
+       JRST    .+4
+       MOVE    D,(A)+1 ;INTERCHANGE THEM
+       MOVEM   D,(B)+1
+       MOVEM   C,(A)+1
+       MOVE    A,B     ;MAKE THE COMPONENT IN "B" THE CURRENT ONE
+       SOJG    E,CLSORT
+       MOVE    E,(P)+1 ;Restore the iteration depth
+       JRST    LLSORT
+\f;ROUTINE TO SORT VECTORS IN NUMERICAL ORDER
+
+VSORT: HLRE    D,VAL1  ;GET COUNT FIELD OF VECTOR
+       IDIV    D,[-2]  ;LENGTH
+       JUMPE   D,GOBACK        ;EMPTY VECTOR?
+       MOVE    E,D     ;SAVE LENGTH IN "E"
+       HRRZ    A,VAL1  ;POINTER TO VECTOR
+       MOVE    B,(A)   ;TYPE OF FIRST COMPONENT
+       CAME    B,$TFIX
+       CAMN    B,$TFLOAT
+       SKIPA
+       JRST    WTYP
+       SOJLE   D,GOBACK        ;IF ONLY ONE COMPONENT THEN FINISHED
+VCOUNT:        ADDI    A,2     ;CHECK NEXT COMPONENT
+       CAME    B,(A)   ;SAME TYPE AS FIRST COMPONENT?
+       JRST    WTYP
+       SOJG    D,VCOUNT        ;CONTINUE WITH NEXT COMPONENT
+
+VVSORT:        SOJE    E,GOBACK        ;FINISHED SORTING?
+       HRRZ    A,VAL1  ;START THIS LOOP OF SORTING AT THE BEGINNING
+       MOVEM   E,(P)+1 ;Save the iteration depth
+CVSORT:        MOVE    C,(A)+3 ;VALUE OF NEXT COMPONENT
+       XCT     O       ;ARE THESE TWO COMPONENTS IN ORDER?
+       JRST    .+4
+       MOVE    D,(A)+1 ;INTERCHANGE THEM
+       MOVEM   D,(A)+3
+       MOVEM   C,(A)+1
+       ADDI    A,2     ;UPDATE THE CURRENT COMPONENT
+       SOJG    E,CVSORT
+       MOVE    E,(P)+1 ;Restore the iteration depth
+       JRST    VVSORT
+
+
+MFUNCTION TIME,SUBR
+       ENTRY 0
+       .RDTIME B,      ;Get time since SYSTEM up
+       MOVSI   A,TFIX
+       JRST    FINIS
+
+
+END
+\f\f\ 3\f
\ No newline at end of file