'*********************************************************************
'* This program solves a Vandermonde linear system:                  *
'*               N      k-1                                          *
'*              S     xi    wi = qk  (k=1,...,N)                     *
'*               i=1                                                 *
'* ----------------------------------------------------------------- *
'* SAMPLE RUN:                                                       *
'*   N=10                                                            *
'*   X(I)=I (I=1,10)                                                 *
'*   Q(I)=1, except Q(9)=30                                          *
'*                                                                   *
'* The solution vector is:                                           *
'*   1    1.004315476190476                                          *
'*   2   -.0381200396825397                                          *
'*   3    .1496031746031746                                          *
'*   4   -.3423611111111111                                          *
'*   5    .5034722222222222                                          *
'*   6   -.4934027777777778                                          *
'*   7    .3222222222222222                                          *
'*   8   -.135218253968254                                           *
'*   9    3.308531746031746D-02                                      *
'*   10   -3.59623015873016D-03                                      *
'*                                                                   *
'* ----------------------------------------------------------------- *
'* Reference: "Numerical Recipes by W.H. Press, B.P. Flannery, S.A.  *
'*             Teukolsky, W.T. Vetterling, Cambridge University      *
'*             Press, 1987"                                          *
'*                                                                   *
'*                             Basic Version By J-P Moreau, Paris    *
'*                                    (www.jpmoreau.fr)              *
'*********************************************************************
'Program Test_Vander
DEFDBL A-H, O-Z
DEFINT I-N

DIM X(10), Q(10), W(10)
  
N = 10
  
DATA 1.0,2.0,3.0,4.0,5.0,6.0,7.0,8.0,9.0,10.0
  
FOR i = 1 TO N
    READ X(i)
NEXT i
    
DATA 1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,30.0,1.0

FOR i = 1 TO N
    READ Q(i)
NEXT i
  
CALL Vander(X(), W(), Q(), N)
  
PRINT
PRINT " The solution vector is:"
  
FOR i = 1 TO N
    PRINT "  "; i; "  "; W(i)
NEXT i

END 'of main program;


SUB Vander (X(), W(), Q(), N)
'---------------------------------------------------------------------
'                                        N     k-1
' Solves the Vandermonde linear system S     xi    wi = qk (k=1,...,N)
'                                        i=1
' Input consists of the vectors X and Q, each of length N, the vector
' W is output.
'---------------------------------------------------------------------
DIM C(N)
'DIM B, ONE, S, T, XX, ZERO AS DOUBLE
ZERO = 0.0
ONE = 1.0
IF N = 1 THEN
    W(1) = Q(1)
ELSE
    FOR i = 1 TO N 'initialize array C
        C(i) = ZERO
    NEXT i
    C(N) = -X(1) 'coefficients of the master polynomial are found by recursion.
    FOR i = 2 TO N
        XX = -X(i)
        FOR j = N + 1 - i TO N - 1
            C(j) = C(j) + XX * C(j + 1)
        NEXT j
        C(N) = C(N) + XX
    NEXT i
    FOR i = 1 TO N 'each subfactor in turn
        XX = X(i)
        T = ONE
        B = ONE
        S = Q(N)
        k = N
        FOR j = 2 TO N 'is synthetically divided,
            k1 = k - 1
            B = C(k) + XX * B
            S = S + Q(k1) * B 'matrix-multiplied by the right-hand side,
            T = XX * T + B
            k = k1
        NEXT j
        W(i) = S / T 'and suppliedwith a denominator.
    NEXT i
END IF
END SUB

'end of file tvander.bas