nxu/fix16.tal

( fix16.tal                                               )
(                                                         )
( use a signed 16-bit short as a fixed point number.      )
(                                                         )
( numbers are interpreted as fractions with an implicit   )
( 256 denominator. the upper byte is signed and           )
( represents the "whole" part of the number, and the      )
( lower byte is unsigned and represents the               )
( "fractional" part of the number.                        )
(                                                         )
( 16-bit fixed point can represent fractional values      )
( in the range -128 <= x < 128. the smallest fraction it  )
( can represent is 1/256, which is about 0.0039.          )
(                                                         )
( SHORT     FRACTION    DECIMAL                           )
( #0000        0/256      0.0000                          )
( #0001        1/256      0.0039                          )
( #0002        2/256      0.0078                          )
( #0040       64/256      0.2500                          )
( #0080      128/256      0.5000                          )
( #0100      256/256      1.0000                          )
( #0700     1792/256      7.0000                          )
( #7f00    32512/256    127.0000                          )
( #7fff    32767/256    127.9961                          )
( #8000   -32768/256   -128.0000                          )
( #8001   -32767/256   -127.9961                          )
( #8100   -32767/256   -127.0000                          )
( #ff00     -256/256     -1.0000                          )
( #ffff       -1/256     -0.0039                          )
(                                                         )
( many 8.8 operations are equivalent to unsigned int16:   )
(   * addition                                            )
(   * subtraction                                         )
( or signed int16:                                        )
(   * comparisons/equality                                )
(                                                         )
( but due to 16-bit truncation multiplication differs...  )
(                                                         )
( x*y = x0*y0 + x0*y1/256 + x1*y0/256 + x1*y1/65536       )
(                                                         )
( since we only have 16-bits:                             )
(   1. we need to drop the 8 high bits from x0*y0         )
(   2. we need to drop the 8 low bits from x1*y1          )
(   3. we need to use all the bits from x0*y1 and x1*y0   )
(                                                         )
( that said, if either x or y is whole (i.e. ends in 00)  )
( then we can just shift that argument right by 8 and use )
( MUL2.                                                   )
(                                                         )
( similarly with division we have:                        )
( )
( x divided by #0100 -> x )
( x divided by #0001 -> x * 256, with likely overfow )
( #abcd divided by #1000 -> #0abc, rounding with d )
( x divided by #0200 -> x >> 1 )
( x divided by #0080 -> x << 1 )

%xyz { ;x16-emit JSR2 #0a #18 DEO }

( useful constants                       )
(                                        )
( to generate your own:                  )
(                                        )
(   1. take true value, e.g. 3.14159...  )
(   2. multiply by 256                   )
(   3. round to nearest whole number     )
(   4. emit hex output                   )
(                                        )
( in python: hex(round(x * 256))         )
%x16-zero      { #0000 } (    0.0        )
%x16-one       { #0100 } (    1.0        )
%x16-two       { #0200 } (    2.0        )
%x16-ten       { #0a00 } (   10.0        )
%x16-hundred   { #6400 } (  100.0        )
%x16-minus-one { #7f00 } (   -1.0        )
%x16-minus-two { #7e00 } (   -2.0        )
%x16-pi/2      { #0192 } (    1.57079... )
%x16-pi        { #0324 } (    3.14159... )
%x16-pi*2      { #0648 } (    6.28318... )
%x16-e         { #02b8 } (    2.71828... )
%x16-phi       { #019e } (    1.61803... )
%x16-sqrt-2    { #016a } (    1.41421... )
%x16-sqrt-3    { #01bb } (    1.73205... )
%x16-epsilon   { #0001 } (    0.00390... )
%x16-minimum   { #8000 } ( -128.0        )
%x16-maximum   { #7fff } (  127.99609... )
%x16-max-whole { #7f00 } (  127.0        )

( useful macros )
%x16-is-non-neg { x16-minimum LTH2 }
%x16-is-neg { x16-maximum GTH2 }

%x16-emit-dec { #30 ADD #18 DEO }

( |0100
    x16-zero xyz
    x16-one xyz
    x16-two xyz
    x16-ten xyz
    x16-hundred xyz
    x16-pi/2 xyz
    x16-pi xyz
    x16-pi*2 xyz
    x16-e xyz
    x16-phi xyz
    x16-sqrt-2 xyz
    x16-sqrt-3 xyz
    x16-epsilon xyz
    #0002 xyz
    #1234 xyz
    #7fff xyz
    #8000 xyz
    #8001 xyz
    #ffff xyz
    #0200 xyz
    #0834 #0000 ;x16-add JSR2 xyz
    #0834 #0834 ;x16-add JSR2 xyz
    #0834 #0834 ;x16-mul JSR2 xyz
    #0834 #0200 ;x16-div JSR2 xyz
    #0100 xyz
    #0003 xyz
    #0100 #0003 ;x16-div JSR2 xyz
    #0834 xyz
    #0080 xyz
    #0834 #0080 ;x16-div JSR2 xyz
    #0834 xyz
    #0300 xyz
    #0834 #0300 ;x16-div JSR2 xyz
    BRK )

@x16-emit ( x* -> )
    DUP2 #8000 EQU2 ,&is-min JCN
    DUP2 #8000 GTH2 ,&is-neg JCN
    SWP DUP #64 LTH ,&<100 JCN
        #64 DIVk DUP x16-emit-dec MUL SUB ,&>=10 JMP
    &is-min POP2
        LIT "- #18 DEO LIT "1 #18 DEO LIT "2 #18 DEO LIT "8 #18 DEO
        LIT ". #18 DEO LIT "0 #18 DEO LIT "0 #18 DEO LIT "0 #18 DEO
        JMP2r
    &is-neg
        LIT "- #18 DEO #ffff EOR2 INC2 ,x16-emit JMP
    &<100 DUP #0a LTH ,&<10 JCN
        &>=10 #0a DIVk DUP x16-emit-dec MUL SUB
    &<10  x16-emit-dec
        LIT '. #18 DEO
    ( emit fractional part )
    #00 SWP ( lo* )
    #000a MUL2 #0100 DIV2k DUP2 NIP x16-emit-dec MUL2 SUB2
    #000a MUL2 #0100 DIV2k DUP2 NIP x16-emit-dec MUL2 SUB2
    #000a MUL2 #0100 DIV2k DUP2 NIP x16-emit-dec MUL2 SUB2
    #000a MUL2 #0100 DIV2k STH2k MUL2 SUB2 #0080 LTH2 ,&no-round JCN INC2r
    &no-round STH2r NIP x16-emit-dec JMP2r

( comparison between x and y. )
(  - ff: x < y )
(  - 00: x = y )
(  - 01: x > y )
@x16-cmp ( x* y* -> c^ )
    STH2k x16-is-neg ,&yn JCN          ( x* [y*] ; ? )
        DUP2 x16-is-non-neg ,&same JCN ( x* [y*] ; y>=0 )
              POP2 POP2r #ff JMP2r     ( -1      ; x<0 y>=0 )
    &yn DUP2 x16-is-neg ,&same JCN     ( x* [y*] ; y<0 )
              POP2 POP2r #01 JMP2r     ( 1       ; x>=0  y<0 )
        &same STH2r ;x16-ucmp JMP2     ( res     ; x<0  y<0 b )

( unsigned comparison between x and y. )
(  - ff: x < y )
(  - 00: x = y )
(  - 01: x > y )
@x16-ucmp ( x* y* -> c^ )
    LTH2k ,&lt JCN GTH2 JMP2r
    &lt POP2 POP2 #ff JMP2r

@x16-eq   ( x* y* ->  x=y^ ) EQU2 JMP2r
@x16-ne   ( x* y* -> x!=0^ ) NEQ2 JMP2r
@x16-lt   ( x* y* ->  x<y^ ) ;x16-cmp JSR2 #ff EQU JMP2r
@x16-lteq ( x* y* ->  x<y^ ) ;x16-cmp JSR2 #01 NEQ JMP2r
@x16-gt   ( x* y* ->  x<y^ ) ;x16-cmp JSR2 #01 EQU JMP2r
@x16-gteq ( x* y* ->  x<y^ ) ;x16-cmp JSR2 #ff NEQ JMP2r

@x16-is-whole ( x* -> bool^ )
    NIP #00 EQU JMP2r

@x16-add ( x* y* -> x+y* )
    ADD2 JMP2r

@x16-sub ( x* y* -> x-y* )
    SUB2 JMP2r

@x16-negate ( x* -> -x* )
    #0000 SWP2 SUB2 JMP2r

@x16-mul ( x* y* -> xy* )
         DUP #00 EQU ,&rhs-whole JCN
    SWP2 DUP #00 EQU ,&rhs-whole JCN
    ,&y3 STR ,&y1 STR ,&x3 STR ,&x1 STR
    LIT2 &x2 00 &x3 00 LIT2 &y2 00 &y3 00 MUL2 #08 SFT2
    LIT2 &x0 00 &x1 00 ,&y2 LDR2          MUL2          ADD2
    ,&x2 LDR2          LIT2 &y0 00 &y1 00 MUL2          ADD2
    ,&x0 LDR2          ,&y0 LDR2 MUL2 #80 SFT2 ADD2 JMP2r
    &rhs-whole #08 SFT2 MUL2 JMP2r

@x16-div ( x* y* -> x/y* )
(    DUP ,&not-whole JCN
        #08 SFT2 DIV2 JMP2r ( since y is whole, x/(y>>8) is correct )
    &not-whole )
    DIV2k STH2k       ( x y x/y  {x/y}         )
    LITr 80 SFT2r     ( x y x/y  {div=(x/y)<<8 )
    OVR2 STH2         ( x y x/y  {y div}       )
    MUL2 SUB2         ( x%y      {y div}       )
    STH2r LIT2r 0100  ( x%y y    {0100 div}    )
    ( we know x%y < y, so start right-shifting y )
    &loop
        DUP2 #0000 EQU2 ,&done JCN
        #01 SFT2 LITr 01 SFT2r  ( rem yi    {shifti div} )
        LTH2k ,&loop JCN        ( rem yi    {shifti div} )
        SWP2 OVR2 SUB2 SWP2     ( rem-yi yi {shifti div} )
        DUP2r ROT2r ADD2r SWP2r ( rem-yi yi {shifti div+shifti} )
        ,&loop JMP              ( rem-yi yi {shifti div+shifti} )
    &done
        POP2 POP2         ( {shiftk div} )
        POP2r STH2r JMP2r ( div )

@x16-quotient ( x* y* -> x//y* )
    DIV2 #80 SFT2 JMP2r

@x16-remainder ( x* y* -> x%y* )
    DIV2k MUL2 SUB2 JMP2r