;
; Calculate and print faculties from 0! to 12!
;
; 2011-10-23 by BlackJack/Civitas
;
        !to "test.prg", cbm

        chrout = $ffd2
        intout = $bdcd
        strout = $ab1e

        i = $02         ; 1 byte
        n = $f7         ; 4 byte

        ; temporary 32 bit numbers
        a = $a5         ; 4 byte
        b = $b0         ; 4 byte
        result = a

        ; Used in print_n routine
        index = $9e     ; 1 byte
        n_str = $57     ; 10 byte

!macro mov32 .source, .dest {
    !for .i, 4 {
        lda .source+.i-1
        sta .dest+.i-1
    }
}

!macro movi8_m32 .val, .dest {
        lda #.val
        !for .i, 4 { sta .dest+.i-1 }
}

!macro rol_n .addr, .n {
    !for .i, .n { rol .addr+.i-1 }
}

!macro rol32 .addr { +rol_n .addr, 4 }

!macro asl32 .addr {
        asl .addr
        +rol_n .addr+1, 3
}

;--------------------------------------
        * = $1000
!zone
start:  ldx #1          ; n = 1
        stx n
        dex
        stx n+1
        stx n+2
        stx n+3
        stx i           ; i = 0

.loop:  lda #0          ; print i
        ldx i
        jsr intout

        lda #<txt       ; print "! = "
        ldy #>txt
        jsr strout

        jsr print_n

        lda #13         ; print newline
        jsr chrout

        ldx i           ; n = n * (i + 1)
        inx
        jsr multiply_n_x
        +mov32 result, n

        inc i           ; i++
        lda i           ; loop finished?
        cmp #13
        bne .loop

        rts

;--------------------------------------
; input: n
; uses: a, index, result, n_str
; calls: divide32_by_10
!zone
print_n
        +mov32 n, a     ; a = n

        lda #0          ; index = 0
        sta index

-       jsr divide32_by_10      ; result, b = divmod(a, 10)

        ldy index       ; n_str[index++] = b
        lda b
        sta n_str,y
        inc index

        lda result      ; result == 0?  -> end of conversion
        !for .i, 3 { ora result+.i }
        beq +

        +mov32 result, a        ; a = result

        jmp -           ; next digit

+       ldy index       ; print digits in n_str backwards
        dey
-       lda n_str,y
        ora #"0"
        jsr chrout
        dey
        bpl -

        rts

;--------------------------------------
; input: a      number to divide by 10
; output: result, b (=remainder)
; uses: b
!zone
divide32_by_10
        +movi8_m32 0, b ; b = 0

        ldx #32         ; all 32 bits
.loop   +asl32 a        ; a = a * 2; b = b * 2 + msb of a
        +rol32 b

        sec             ; mod1[1]/mod2[1]/Y/A = b - 10
        lda b
        sbc #10
        sta .mod1+1
        lda b+1
        sbc #0
        sta .mod2+1
        lda b+2
        sbc #0
        tay
        lda b+3
        sbc #0
        bcc +           ; if b<10 don't store b-10 in b, don't set lsb in result

        sta b+3         ; b = b - 10
        sty b+2
.mod2:  lda #0
        sta b+1
.mod1:  lda #0
        sta b
        inc result

+       dex             ; next bit
        bne .loop
        rts

;--------------------------------------
; input: n, X
; output: result
; uses: b
!zone
multiply_n_x
        +movi8_m32 0, result    ; result = 0
        +mov32 n, b             ; b = n

        txa
        ldx #8
.loop   lsr
        bcc +           ; if bit not set skip add

        tay             ; result += b
        clc
        !for .i, 4 {
            lda result+.i-1
            adc b+.i-1
            sta result+.i-1
        }
        tya

+       +asl32 b        ; b *= 2

        dex             ; next bit
        bne .loop

        rts

;--------------------------------------
txt:    !pet "! = ", 0