Copyright (C) 2020 Andreas Kloeckner

In [5]:

```
import numpy as np
import numpy.linalg as la
```

We start by choosing our *quadrature nodes*, the maximum degree which will be exact, as well as the interval $(a,b)$ on which we integrate:

In [6]:

```
nodes = [0, 1]
#nodes = [0, 0.5, 1]
#nodes = [3, 3.5, 4]
#nodes = [0, 1, 2]
#nodes = np.linspace(0, 1, 12)
max_degree = len(nodes)-1
a = nodes[0]
b = nodes[-1]
```

Next, we compute the transpose of the Vandermonde matrix $V^T$ and the integrals $\int_a^b x^i$ as `rhs`

:

In [13]:

```
nodes = np.array(nodes)
powers = np.arange(max_degree+1)
Vt = nodes ** powers.reshape(-1, 1)
rhs = 1/(powers+1) * (b**(powers+1) - a**(powers+1))
if len(nodes) <= 4:
print(Vt)
```

[[1 1] [0 1]]

Set up the linear system for the weights:

$$ \begin{align*} \alpha_0 x_0^0 + \cdots + \alpha_{n-1} x_{n-1}^{0} &= \int_a^b x^0\\ \vdots &= \vdots \\ \alpha_0 x_0^{n-1} + \cdots + \alpha_{n-1} x_{n-1}^{n-1} &= \int_a^b x^{n-1} \end{align*} $$In [9]:

```
weights = la.solve(Vt, rhs)
print(weights)
```

[ 0.5 0.5]

Now we test our quadrature rule by integrating the monomials $\int_a^b x^i dx$ and comparing quadrature results to the true answers:

In [12]:

```
for i in range(len(nodes) + 1):
approx = weights @ nodes**i
true = 1/(i+1)*(b**(i+1) - a**(i+1))
print("Error at degree %d: %g" % (i, approx-true))
```

Error at degree 0: 0 Error at degree 1: 0 Error at degree 2: 0.166667

In [24]:

```
```