{ "cells": [ { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Energy System Modelling - Solutions to Tutorial II" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "skip" } }, "source": [ "**Imports**" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "slideshow": { "slide_type": "skip" } }, "outputs": [], "source": [ "import numpy as np\n", "import numpy.linalg\n", "import pandas as pd\n", "import matplotlib.pyplot as plt\n", "%matplotlib inline" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "fragment" } }, "source": [ "***\n", "## Problem II.1" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "***\n", "**(a) Compile the nodes list and the edge list.**\n", "> **Remark:** While graph-theoretically both lists are unordered sets, let's agree on an ordering now which can serve as basis for the matrices in the following exercises: we sort everything in ascending numerical order, i.e.\\ node 1 before node 2 and edge (1,2) before (1,4) before (2,3)." ] }, { "cell_type": "code", "execution_count": 142, "metadata": { "slideshow": { "slide_type": "fragment" } }, "outputs": [], "source": [ "nodes = [0, 1, 2, 3, 4, 5]\n", "edges = [(0, 1), (1, 2), (1, 4),\n", " (1, 5), (2, 3), (2, 4),\n", " (3, 4), (4, 5)]" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "***\n", "**(b) Determine the order and the size of the network.**" ] }, { "cell_type": "code", "execution_count": 143, "metadata": { "slideshow": { "slide_type": "fragment" } }, "outputs": [], "source": [ "N = len(nodes)\n", "E = len(edges)" ] }, { "cell_type": "code", "execution_count": 144, "metadata": { "slideshow": { "slide_type": "fragment" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Order: 6\n", "Size: 8\n" ] } ], "source": [ "print(\"Order: {}\\nSize: {}\".format(N, E))" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "***\n", "**(c) Compute the adjacency matrix $A$ and check that it is symmetric.**" ] }, { "cell_type": "code", "execution_count": 145, "metadata": { "slideshow": { "slide_type": "fragment" } }, "outputs": [ { "data": { "text/plain": [ "array([[0., 1., 0., 0., 0., 0.],\n", " [1., 0., 1., 0., 1., 1.],\n", " [0., 1., 0., 1., 1., 0.],\n", " [0., 0., 1., 0., 1., 0.],\n", " [0., 1., 1., 1., 0., 1.],\n", " [0., 1., 0., 0., 1., 0.]])" ] }, "execution_count": 145, "metadata": {}, "output_type": "execute_result" } ], "source": [ "A = np.zeros((N, N))\n", "\n", "for u, v in edges:\n", " A[u, v] += 1\n", " A[v, u] += 1\n", "\n", "A" ] }, { "cell_type": "code", "execution_count": 146, "metadata": { "slideshow": { "slide_type": "fragment" } }, "outputs": [ { "data": { "text/plain": [ "True" ] }, "execution_count": 146, "metadata": {}, "output_type": "execute_result" } ], "source": [ "(A == A.T).all()" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "***\n", "**(d) Find the $k_n$ of each node $n$ and compute the average degree of the network.**" ] }, { "cell_type": "code", "execution_count": 147, "metadata": { "slideshow": { "slide_type": "fragment" } }, "outputs": [ { "data": { "text/plain": [ "array([1., 4., 3., 2., 4., 2.])" ] }, "execution_count": 147, "metadata": {}, "output_type": "execute_result" } ], "source": [ "k = A.sum(axis=1)\n", "k" ] }, { "cell_type": "code", "execution_count": 148, "metadata": { "slideshow": { "slide_type": "fragment" } }, "outputs": [ { "data": { "text/plain": [ "2.6666666666666665" ] }, "execution_count": 148, "metadata": {}, "output_type": "execute_result" } ], "source": [ "k.mean()" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "***\n", "**(e) Determine the incidence matrix $K$ by assuming the links are always directed from smaller-numbered node to larger-numbered node, i.e.\\ from node 2 to node 3, instead of from 3 to 2.**" ] }, { "cell_type": "code", "execution_count": 150, "metadata": { "slideshow": { "slide_type": "fragment" } }, "outputs": [ { "data": { "text/plain": [ "array([[ 1., 0., 0., 0., 0., 0., 0., 0.],\n", " [-1., 1., 1., 1., 0., 0., 0., 0.],\n", " [ 0., -1., 0., 0., 1., 1., 0., 0.],\n", " [ 0., 0., 0., 0., -1., 0., 1., 0.],\n", " [ 0., 0., -1., 0., 0., -1., -1., 1.],\n", " [ 0., 0., 0., -1., 0., 0., 0., -1.]])" ] }, "execution_count": 150, "metadata": {}, "output_type": "execute_result" } ], "source": [ "K = np.zeros((N,E))\n", "\n", "for i, (u, v) in enumerate(edges):\n", " K[u,i] = 1\n", " K[v,i] = -1\n", " \n", "K" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "***\n", "**(f) Compute the Laplacian $L$ of the network using $k_n$ and $A$. Remember that the Laplacian can also be computed as $L=KK^T$ and check that the two definitions agree.**" ] }, { "cell_type": "code", "execution_count": 151, "metadata": { "slideshow": { "slide_type": "fragment" } }, "outputs": [], "source": [ "D = np.diag(A.sum(axis=1))" ] }, { "cell_type": "code", "execution_count": 152, "metadata": { "slideshow": { "slide_type": "fragment" } }, "outputs": [ { "data": { "text/plain": [ "array([[ 1., -1., 0., 0., 0., 0.],\n", " [-1., 4., -1., 0., -1., -1.],\n", " [ 0., -1., 3., -1., -1., 0.],\n", " [ 0., 0., -1., 2., -1., 0.],\n", " [ 0., -1., -1., -1., 4., -1.],\n", " [ 0., -1., 0., 0., -1., 2.]])" ] }, "execution_count": 152, "metadata": {}, "output_type": "execute_result" } ], "source": [ "L = D - A\n", "L" ] }, { "cell_type": "code", "execution_count": 153, "metadata": { "slideshow": { "slide_type": "fragment" } }, "outputs": [ { "data": { "text/plain": [ "True" ] }, "execution_count": 153, "metadata": {}, "output_type": "execute_result" } ], "source": [ "np.array_equal(L, K.dot(K.T))" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "***\n", "**(g) Find the diameter of the network by looking at the graph.**" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "fragment" } }, "source": [ "By looking: Between nodes $0$ and $3$, f.ex. $0 \\to 1 \\to 2 \\to 3$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "***\n", "## Problem II.3\n", "\n", "If you map the nodes to `0=DK, 1=DE, 2=CH, 3=IT, 4=AT,5=CZ` the network represents a small part of the European electricity network (albeit very simplified). On the [course homepage](https://nworbmot.org/courses/complex_renewable_energy_networks/), you can find the power imbalance time series for the six countries for January 2017 in hourly MW in the file `imbalance.csv`. They have been derived from physical flows as published by [ENTSO-E](https://transparency.entsoe.eu/transmission-domain/physicalFlow/show)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "The linear power flow is given by\n", "$$p_i = \\sum_j \\tilde{L}_{i,j}\\theta_j \\qquad \\text{and} \\qquad f_l = \\frac{1}{x_l} \\sum_i K_{i,l}\\theta_i, \\qquad \\text{where} \\qquad \\tilde{L}_{i,j}= \\sum_l = K_{i,l}\\frac{1}{x_l} K_{j,l}$$\n", "is the weighted Laplacian. For simplicity, we assume identity reactance on all links $x_l = 1$." ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "***\n", "**Read data**" ] }, { "cell_type": "code", "execution_count": 154, "metadata": { "slideshow": { "slide_type": "fragment" } }, "outputs": [ { "data": { "text/html": [ "
\n", " | DK | \n", "DE | \n", "CH | \n", "IT | \n", "AT | \n", "CZ | \n", "
---|---|---|---|---|---|---|
DateTime | \n", "\n", " | \n", " | \n", " | \n", " | \n", " | \n", " |
2017-01-01 00:00:00 | \n", "-79.73 | \n", "3486.25 | \n", "-3403.37 | \n", "-97.5 | \n", "-2867.33 | \n", "2961.68 | \n", "
2017-01-01 01:00:00 | \n", "-109.14 | \n", "4248.17 | \n", "-3829.59 | \n", "-105.9 | \n", "-3149.98 | \n", "2946.44 | \n", "
2017-01-01 02:00:00 | \n", "-101.49 | \n", "5314.56 | \n", "-4131.66 | \n", "-199.2 | \n", "-3646.87 | \n", "2764.66 | \n", "
2017-01-01 03:00:00 | \n", "-197.96 | \n", "5450.71 | \n", "-4009.61 | \n", "-656.0 | \n", "-3345.56 | \n", "2758.42 | \n", "
2017-01-01 04:00:00 | \n", "-526.15 | \n", "5695.72 | \n", "-3924.43 | \n", "-646.2 | \n", "-3393.12 | \n", "2794.18 | \n", "