In [ ]:
import numpy
%matplotlib inline
import matplotlib.pyplot

Excersises

1.

Implemet softmax. Softmax is a vector-vector function such that $$ x_i \mapsto \frac{\exp(x_i)}{\sum_{j=1}^n \exp(x_j)} $$ Avoid using for loops, use vectorization. The solution below is a bad solution.


In [ ]:
def softmax(X):
    X = numpy.array(X)
    Y = numpy.exp(X)
    return Y/Y.sum()
print(softmax([-1,0,1]))

2.

Write a function which has one parameter, a 2D array and it returns a vector of row-wise Euclidean norms of the input. Use numpy operations and vectorization, avoid for loops. The solution below is a bad solution.


In [ ]:
def rowwise_norm(X):
    XX = numpy.array(X)**2
    return numpy.sqrt(XX.sum(1))

X = numpy.arange(5)[:, None]*numpy.ones((5, 3));
print(X)
print(rowwise_norm(X))
print(rowwise_norm([[1], [-1], [1], [-1]]))

3.

Write a function which has one parameter, a positive integer $n$, and returns an $n\times n$ array of $\pm1$ values like a chessboard: $M_{i,j} = (-1)^{i+j}$.

Use numpy operations, avoid for loops!

First just a vector of laternating minus ones and ones.


In [ ]:
print((numpy.arange(5) % 2) * 2 - 1)
print("or")
print((-1)**numpy.arange(5))
print("rather")
print(-(-1)**numpy.arange(5))

In [ ]:
def chessboard(n):
    a = -(-1)**numpy.arange(n)
    return a.reshape((-1, 1)) * a.reshape((1, -1))
    return a[:, None] * a[None, :] # the same like above, just shorter

print(chessboard(5))

4.*

Write a function which numerically derivates a $\mathbb{R}\mapsto\mathbb{R}$ function. Use the forward finite difference.

The input is a 1D array of function values, and optionally a 1D vector of abscissa values. If not provided then the abscissa values are unit steps.

The result is a 1D array with the length of one less than the input array.

Use numpy operations instead of for loop in contrast to the solution below.


In [ ]:
def derivate(f, x=None):
    if x is None:
        x = numpy.arange(len(f))

    df = f[:-1] - f[1:] # vector of differences of y values
    dx = x[:-1] - x[1:] # vector of differences of x (abscissa) values
    return df/dx

print(derivate(numpy.arange(10)**2))

In [ ]:
x=numpy.arange(0,1,0.1)
matplotlib.pyplot.plot(x, x**2)
matplotlib.pyplot.plot(x[:-1], derivate(x**2, x))

4. b)

Make a random (1D) grid and use that as abscissa values!


In [ ]:
x = numpy.sort(numpy.random.rand(100))
x = numpy.concatenate([[0], x, [1]]) # the grid contains 0, a hundered random numbers, and 1 in an increasing order
f = x**2
matplotlib.pyplot.plot(x, f)
matplotlib.pyplot.plot(x[:-1], derivate(f, x))

5.*

Implement the Horner's method for evaluating polynomials. The first input is a 1D array of numbers, the coefficients, from the constant coefficient to the highest order coefficent. The second input is the variable $x$ to subsitute. The function should work for all type of variables: numbers, arrays; the output should be the same type array as the input, containing the elementwise polynomial values.


In [ ]:
def horner(C, x):
    y = numpy.zeros_like(x)
    for c in C[::-1]:
        y = y*x + c
    return y

In [ ]:
C = [2, 0, 1] # 2 + 0*x + 1*x^2
print(horner(C, 3)) # 2 + 3^2
print(horner(C, [3, 3]))
print(horner(C, numpy.arange(9).reshape((3,3))))

With a slight modofication, you can implement matrix polinomials!


In [ ]:
def matrix_horner(C, x):
    assert(len(x.shape) == 2) # matrices
    assert(x.shape[0] == x.shape[1]) # square matrices
    eye = numpy.eye(x.shape[0], dtype=x.dtype)
    
    y = numpy.zeros_like(x)
    for c in C[::-1]:
        y = y.dot(x) + c*eye
    return y

In [ ]:
C = [2, 0, 1] # 2 + 0*x + 1*x^2
M = numpy.arange(9).reshape((3,3))
print(M)
print()
print(M.dot(M))
print()
print(2*numpy.eye(3) + M.dot(M))
print(matrix_horner(C, M))

6.*

Plot the $z\mapsto \exp(z)$ complex function on $[-2, 2]\times i [-2, 2]$. Use matplotlib.pyplot.imshow and the red and green color channgels for real and imaginary parts.


In [ ]:
z = -1j*numpy.linspace(-2,2,100)[:, None] + numpy.linspace(-2,2,100)[None, :] # complex grid
expz = numpy.exp(z) # function values

real_part = expz.real
imag_part = expz.imag

# for scaling to [0,1]
max_real = real_part.max()
min_real = real_part.min()
max_imag = imag_part.max()
min_imag = imag_part.min()

# RGB
colors = numpy.stack([
                        (real_part - min_real) / (max_real - min_real),
                        (imag_part - min_imag) / (max_imag - min_imag),
                        numpy.zeros_like(real_part)
                    ], axis=2)
# plotting
matplotlib.pyplot.imshow(colors)

In [ ]: