using Plots, LinearAlgebra

Vandermonde Matrix

The Vandermonde matrix naturally arises in approximating a function using a polynomial:

$$ f(x) = a_0 + a_1 x + a_2 x^2 + \ldots + a_{N-1} x^{N-1}. $$

By sampling the functions at $N$ distinct points, the coefficients $a_n$'s satisfy

$$ \begin{aligned} V \begin{bmatrix} a_0 \\ a_1 \\ \ldots \\ a_{N-1} \end{bmatrix} = \begin{bmatrix} f(x_0) \\ f(x_1) \\ \ldots \\ f(x_{N-1}) \end{bmatrix} \end{aligned}, $$

where

$$ \begin{aligned} V = \begin{bmatrix} 1 & x_0 & x_0^2 & \ldots & x_0^{N-1} \\ 1 & x_1 & x_1^2 & \ldots & x_1^{N-1} \\ 1 & x_2 & x_2^2 & \ldots & x_2^{N-1} \\ & & \ldots & & \\ 1 & x_{N-1} & x_{N-1}^2 & \ldots & x_{N-1}^{N-1} \end{bmatrix} \end{aligned} $$

function Vandermonde(xv)
    
    ndim = length(xv)
    vmat = zeros(ndim, ndim)
    
    for i1 in 1:ndim, i2 in 1:ndim
			vmat[i1,i2] = xv[i1]^(i2-1)
	end
    return vmat
    
end

function Vinterp(f, xval; a=-1, b=1, ndim=25)

    xv = range(a, b, length=ndim)
    coeffs = Vandermonde(xv) \ f.(xv)
    ret = sum( [coeffs[j] * xval^(j-1) for j in 1:ndim] )
    
    return ret
    
end

Vinterp (generic function with 1 method)

f = x -> 1/(1+x^2)

_f1 = x -> Vinterp(f, x, ndim=5)
xv = range(-1.2, 1.2, length=20)

p1 = plot(xv, f, label="exact")
scatter!(xv, _f1, label="interpolation (5th)")

Quadrature Formula

$$ \int_{-1}^{1} dx f(x) \approx \sum_{n={\rm even}} \frac{2}{n+1} a_n = \sum_{m} y_m \sum_{n={\rm even}} \frac{2}{n+1} V^{-1}_{nm} $$

vv = Vandermonde(range(-1, 1, length=3))

function quad_coeff(i1; norder=3)

    vv = Vandermonde(range(-1, 1, length=norder))
    
    ret = 0.0

    for j1 in 1:3
        if j1%2 == 0
            continue
        else
            ret += inv(vv)[j1, i1] * (2/(j1))
        end
    end

    return ret
end

quad_coeff.(1:3, norder=3)