PROGRAM D2r10

! Driver for sparse

LIBRARY "sparse", "asub"

DECLARE PUBLIC n                  ! COMMON with asub and atsub, below

LET n = 20
DIM b(0), x(0), bcmp(0)
MAT redim b(n), x(n), bcmp(n)

CLEAR
LET m = n
MAT x = zer
MAT b = con
LET b(1) = 3
LET b(n) = -1

CALL sparse (b(), n, dum, dum, x(), rsq)

PRINT "Sum-squared residual:",
PRINT using "#.######^^^^": rsq
PRINT
PRINT "Solution vector:"
FOR i = 1 to n
    PRINT using "----#.######": x(i);
    IF mod(i, 5) = 0 then PRINT
NEXT i
PRINT
PRINT

CALL asub (x(), bcmp())

PRINT "Press Enter to continue ..."
LINE INPUT dum$

PRINT "Test of solution vector:"
PRINT "      a*x", "      b"
FOR i = 1 to n
    PRINT using "----#.######": bcmp(i), b(i)
NEXT i
PRINT

END

MODULE sparsesubs

    PUBLIC n                      ! Set in the main program.

    SUB asub (xin(), xout())

        LET xout(1) = xin(1) + 2 * xin(2)
        LET xout(n) = -2 * xin(n - 1) + xin(n)
        FOR i = 2 to n - 1
            LET xout(i) = -2 * xin(i - 1) + xin(i) + 2 * xin(i + 1)
        NEXT i

    END SUB

    SUB atsub (xin(), xout())

        LET xout(1) = xin(1) - 2 * xin(2)
        LET xout(n) = 2 * xin(n - 1) + xin(n)
        FOR i = 2 to n - 1
            LET xout(i) = 2 * xin(i - 1) + xin(i) - 2 * xin(i + 1)
        NEXT i

    END SUB

END MODULE