complete exercise 2

exercise2/exercise2.typ

@@ -0,0 +1,404 @@
+#import "@preview/cetz:0.3.2";
+#import "@preview/cetz-plot:0.1.1": plot
+#import "@preview/physica:0.9.4": *
+#import "@preview/plotsy-3d:0.1.0": plot-3d-parametric-surface
+#import "@preview/fletcher:0.5.4" as fletcher: diagram, edge, node
+
+#set page(paper: "a4", margin: (x: 2.6cm, y: 2.8cm), numbering: "1 : 1")
+#set par(justify: true, leading: 0.52em)
+
+#let FONT_SIZE = 18pt;
+#set text(font: "FreeSerif", size: FONT_SIZE, lang: "us")
+#show math.equation: set text(font: "Euler Math", size: (FONT_SIZE * 1.0), lang: "en")
+
+#set heading(numbering: none)
+#show heading.where(level: 1): it => {
+  rect(inset: FONT_SIZE / 2)[#it]
+}
+// #show heading.where(level: 2): it => [
+//   #set align(left)
+//   #set text(size: FONT_SIZE * 1.1, weight: "semibold")
+//   #(it.body)
+// ]
+
+#align(center)[
+  #text(size: FONT_SIZE * 2, weight: "bold")[#underline[exercise 1]]
+]
+
+these are my solutions to the second exercise set of TMA4135.
+
+there should be a python source code file attached to this deliverable.
+this file can be used to generate the images shown in problem 2.
+
+this document was created using
+#link("https://typst.app/")[#text(blue.darken(5%))[typst]].
+
+#v(42pt)
+
+#outline(title: none)
+
+#v(42pt)
+
+
+= problem 1
+
+define the _cardinal functions_ as
+$
+  ell_i(x) & := product_(j=0\ j!=i)^n (x - x_j) / (x_i - x_j) \
+           & = (x - x_0)/(x_i - x_0) dot dots.c dot
+             (x - x_(i-1))/(x_i - x_(i-1)) dot (x - x_(i+1))/(x_i - x_(i+1))
+             dot dots.c dot (x - x_n)/(x_i - x_n),
+$
+for $i = 0, ..., n$.
+
+
+== a)
+
+given $x_0 = -3, x_1 = -1, #[and] x_2 = 5$, then
+#table(
+  $
+    l_0 & = product_(j=0\ j!=i)^n (x - x_j)/(x_0 - x_j) \
+        & = (x - x_1)/(x_0 - x_1) dot (x - x_2)/(x_0 - x_2) \
+        & = (x + 1)/(-2) dot (x - 5)/(-8) \
+        & = ((x + 1)(x - 5))/16 \
+        & = (x^2 - 4x - 5)/16
+  $,
+  $
+    l_1 & = product_(j=0\ j!=i)^n (x - x_j)/(x_1 - x_j) \
+        & = (x - x_0)/(x_1 - x_0) dot (x - x_2)/(x_1 - x_2) \
+        & = (x + 3)/(2) dot (x - 5)/(-6) \
+        & = ((x + 3)(x - 5))/(-12) \
+        & = (-x^2 + 2x + 15)/12
+  $,
+  $
+    l_2 & = product_(j=0\ j!=i)^n (x - x_j)/(x_2 - x_j) \
+        & = (x - x_0)/(x_2 - x_0) dot (x - x_1)/(x_2 - x_1) \
+        & = (x + 3)/(8) dot (x + 1)/(6) \
+        & = ((x + 3)(x + 1))/48 \
+        & = (x^2 + 4x + 3)/48
+  $,
+
+  column-gutter: 20%,
+  stroke: none,
+  columns: 3,
+)
+
+
+== b)
+
+we can easily verify these cardinal functions
+$
+  l_0(x_0) & = ((-3)^2 - 4(-3) - 5)/16 = (9 + 12 - 5)/16 = 1 \
+  l_0(x_1) & = ((-1)^2 - 4(-1) - 5)/16 = (1 + 4 - 5)/16 = 0 \
+  l_0(x_2) & = (5^2 - 4 dot 5 - 5)/16 = (25 - 20 - 5)/16 = 0 \
+           \
+           \
+  l_1(x_0) & = (-(-3)^2 + 2(-3) + 15)/12 = (-9 - 6 + 15)/12 = 0 \
+  l_1(x_1) & = (-(-1)^2 + 2(-1) + 15)/12 = (-1 -2 + 15)/12 = 1 \
+  l_1(x_2) & = (-5^2 + 2 dot 5 + 15)/12 = (-25 + 10 + 15)/12 = 0 \
+           \
+           \
+  l_2(x_0) & = ((-3)^2 + 4(-3) + 3)/48 = (9 - 12 + 3)/48 = 0 \
+  l_2(x_1) & = ((-1)^2 + 4(-1) + 3)/48 = (1 -4 + 3)/48 = 0 \
+  l_2(x_2) & = (5^2 + 4 dot 5 + 3)/48 = (25 + 20 + 3)/48 = 1
+$
+
+thus we see that they follow the desired pattern.
+
+
+== c)
+
+assume we have $n + 1$ data points, then
+$
+  p_n(x) := sum_(i=0)^n y_i ell_i(x)
+$
+must be the (unique) interpolation polynomial for points $(x_i, y_i)$.
+
+to show this, we first must convince ourselves that $deg(ell_i (x)) = n$.
+recall the definition
+$
+     ell_i(x) & := product_(j=0\ j!=i)^n (x - x_j) / (x_i - x_j) \
+  => ell_i(x) & = product_(j=0\ j!=i)^n f(x)
+$
+where
+$
+  f(x) = (x - x_j) / (x_i - x_j) \
+  => deg(f(x)) = 1
+$
+
+the product contains exactly $n$ factors, because we iterate over the $n + 1$
+data points, minus the one for $i = j$. thus, $deg(ell_i(x)) = n$.
+
+lastly, we must verify that the interpolation condition holds for $p_n(x)$
+$
+  p_n(x_i) = y_i, quad i = 0, ..., n.
+$
+
+we can expand the definition to show this
+$
+  p_n(x_i) & = sum_(j=0)^n y_j ell_j(x_i) \
+           & = sum_(j=0)^n y_j product_(k=0\ k!=j)^n (x_i - x_k) / (x_j - x_k) \
+           & = sum_(j=0)^n y_j dot g(i, j)
+$
+where
+$
+  g(i, j) := cases(
+    1 & quad "if" i = j,
+    0 & quad "otherwise"
+  )
+$
+
+thus we get
+$
+  p_n(x_i) & = (sum_(j=0\ j!=i)^n y_j dot 0) + (sum_(j=0\ j=i)^n y_j dot 1) \
+           & = (0) + (y_i) = y_i
+$
+
+as such, we can see that $p_n(x)$ is the interpolation polynomial of the given
+data points.
+
+
+== d)
+
+recall from a) that given $x_0 = -3, x_1 = -1, #[and] x_2 = 5$, then
+#table(
+  $
+    l_0 = (x^2 - 4x - 5)/16
+  $,
+  $
+    l_1 = (-x^2 + 2x + 15)/12
+  $,
+  $
+    l_2 = (x^2 + 4x + 3)/48
+  $,
+
+  column-gutter: 20%,
+  stroke: none,
+  columns: 3,
+)
+
+let $y_0 = 8, y_1 = -4, #[and] y_2 = 12$. then using our proven formula from c),
+we can compute the interpolation polynomial
+$
+  p(x) & = sum_(j=0)^2 y_j ell_j(x) \
+       & = y_0 ell_0(x) + y_1 ell_1(x) + y_2 ell_2(x) \
+       & = (x^2 - 4x - 5)/2 + (x^2 - 2x - 15)/3 + (x^2 + 4x + 3)/4 \
+       & = (6x^2 - 24x - 30 + 4x^2 - 8x - 60 + 3x^2 + 12x + 9)/12 \
+       & = (13x^2 - 20x - 81)/12.
+$
+
+
+= problem 2
+
+== a)
+
+#image("x.png", width: 80%)
+
+== b)
+
+#image("y.png", width: 80%)
+
+== c) & d)
+
+#image("path.png", width: 80%)
+
+#pagebreak()
+
+= problem 3
+
+let $y = f(x) := 2 + sin(x) cos(x)$ with
+#[
+  #set align(center)
+  #set table(
+    stroke: (x, y) => {
+      if y == 0 { (bottom: luma(180)) }
+      if x == 0 { (right: luma(180)) }
+    },
+    inset: 10% + 5pt,
+  )
+  #table(
+    $i$, $0$, $1$, $2$,
+    $x_i$, $0$, $pi/6$, $pi/3$,
+    $y_i$, $2$, $2 + sqrt(3)/4$, $2 + sqrt(3)/4$,
+    columns: 4,
+  )
+]
+
+and define
+$
+  p(x) := a_1 + a_2 sin x + a_3 sin 2x
+$
+where $p(x_i) = y_i #[for] i in {0, 1, 2} #[and] a_1, a_2, a_3 in RR$.
+
+== a)
+
+from some quick analysis, we can see that
+$
+  a_1 = y_1 = 2 => \
+  sqrt(3)/4 = a_2/2 + sqrt(3) a_3 / 2 and \
+  sqrt(3)/4 = sqrt(3) a_2 / 2 + sqrt(3) a_3 / 2 => \
+  a_2/2 + sqrt(3) a_3 / 2 = sqrt(3) a_2 / 2 + sqrt(3) a_3 / 2 \
+  a_2 = sqrt(3) a_2 => a_2 = 0 => \
+  a_3 = 1/2
+$
+
+so we get $p(x) = 2 + 1/2 sin 2x$.
+
+== b)
+
+now, we can show that $p(x) = f(x)$, i.e. the error $e(x) = 0$
+
+recall that
+$
+  sin 2x = 2 sin x cos x
+$
+
+we then get
+$
+  p(x) & = 2 + 1/2 sin 2x \
+       & = 2 + 1/2 (2 sin x cos x) \
+       & = 2 + sin x cos x = f(x)
+$
+as expected.
+
+= problem 4
+
+define
+$
+  c_0 := y_0, quad
+  p_0(x) := c_0
+$
+and recursively
+$
+      c_k & := (y_k - p_(k-1)(x_k))/(product_(i=0)^(k-1) (x_k - x_i)), \
+  p_k (x) & := p_(k-1)(x) + c_k product_(i=0)^(k-1) (x-x_i),
+$
+for $k = 1, ..., n$.
+
+== a)
+
+given data points
+#[
+  #set align(center)
+  #set table(
+    stroke: (x, y) => {
+      if y == 0 { (bottom: luma(180)) }
+      if x == 0 { (right: luma(180)) }
+    },
+    inset: 10% + 5pt,
+  )
+  #table(
+    $i$, $0$, $1$, $2$,
+    $x_i$, $-1$, $1$, $5$,
+    $y_i$, $6$, $-2$, $4$,
+    columns: 4,
+  )
+]
+
+compute $p_2(x)$.
+
+first we need
+$
+  p_1(x) & = p_0(x) + (y_1 - p_0(x_1))/(x_1 - x_0) dot (x - x_0) \
+         & = 6 + (-2 - 6)/(1 - (-1)) dot (x - (-1)) \
+         & = 6 - 4 (x + 1) = -4x + 2
+$
+then use $p_1(x)$ in
+$
+  p_2(x) & = p_1(x) + (y_2 - p_1(x_2))/((x_2 - x_0)(x_2 - x_1))
+           dot (x - x_0)(x - x_1) \
+         & = -4x + 2 + (4 - (-4)dot 4 + 2)/((5-(-1))(5-1))
+           dot (x - (-1))(x - 1) \
+         & = -4x + 2 + 11/12 dot (x + 1)(x - 1) \
+         & = 11/12 x^2 - 4x + 13/12
+$
+
+#linebreak()
+
+== b)
+
+$p_k (x)$ is the interpolation polynomial for the data points, because
+#[ #set math.equation(numbering: "(1)", supplement: none)
+  $
+    deg(p_k (x)) = k quad forall quad k in NN quad
+  $ <p1>
+  and
+  $
+    p_k (x_i) = y_i quad forall quad i = 0, ..., k quad
+  $ <p2>
+]
+
+to show @p1, we can see from the definition that the product has degree $k$,
+because there are $k$ factors containing $x$. by inductive argument, we can
+trivially see that $p_(k-1) (x)$ must have a degree equal to $k-1$. as such,
+$
+  deg(p_k (x)) = k
+$
+
+next, (@p2) can be shown by a splitting the problem into two cases, one of them
+requiring induction, otherwise just direct proofs. observe
+$
+  p_k (x_k) & = p_(k-1) (x_k) + c_k product_(i=0)^(k-1) (x_k - x_i) \
+            & = cancel(p_(k-1) (x_k)) + (y_k - cancel(p_(k-1) (x_k)))/cancel(
+                product_(i=0)^(k-1) (x_k
+                  - x_i)
+              ) dot cancel(product_(i=0)^(k-1) (x_k - x_i)) \
+            & = y_k.
+$
+thus (@p2) holds for $x_k$. to show that it holds for $x_j$ where $j < k$
+requires that we analyze the product in the definition:
+$
+  product_(i=0)^(k-1) (x_j - x_i)
+$
+but since $j in {0, ..., k - 1}$, the product must contain a $0$-factor for any
+chosen $j$, thus canceling the entire second term of the definition, leaving us
+with
+$
+  p_k (x_j) & = p_(k-1) (x_j).
+$
+then we need to show that
+$
+  p_(k-1) (x_j) & = y_j quad forall quad j in {0, ..., k - 1}
+$
+
+let us perform an inductive argument. the hypothesis is
+$
+  p_0 (x_j) = y_j
+$
+
+this is wrong.
+
+== c)
+
+we can use the points from a) to form a scheme
+
+#[
+  #set align(center)
+  #set table(
+    stroke: (x, y) => {
+      if x == 0 { (right: luma(180)) }
+    },
+    inset: 10% + 5pt,
+  )
+  #table(
+    $x_0 = -1$, $y_(0,0) = 6$, none, none,
+    none, none, $y_(0,1) = -4$, none,
+    $x_1 = 1$, $y_(1,1) = -2$, none, $y_(0, 2) = 13 slash 12$,
+    none, none, $y_(1,2) = 3 slash 2$, none,
+    $x_2 = 5$, $y_(2,2) = 4$, none, none,
+    columns: 4,
+  )
+]
+
+then we can form the interpolation polynomial according to the formula
+$
+  p_k (x) = sum_(i=0)^k y_(0, i) product_(j=0)^(i-1) (x-x_j)
+$
+such that for $n + 1 = 3$ we get
+$
+  p_n (x) = p_2(x) & = 6 - 4(x - (-1)) + 13/12 (x - (-1))(x - 1) \
+                   & = 6 - 4x - 4 + 13/12 (x^2 - 1) \
+                   & = 13/12 x^2 - 4x + 11/12
+$
+which is _almost_ what we got in a) suggesting that i've made a slight error in
+either task.

exercise2/main.py

@@ -0,0 +1,55 @@
+from sympy import prod, cos, sin
+from sympy.abc import t
+from sympy.plotting import plot, plot_parametric
+
+t_i = [0, 0.8976, 1.7952, 2.6928, 3.5904, 4.4880, 5.3856, 6.2832]
+x_i = [1, 1.5984, -0.6564, -1.6828, -0.1191, 0.2114, -0.3514, 1]
+y_i = [0, 0.7818, 0.9750, 0.4339, -0.4339, -0.975, -0.7818, 0]
+n = len(x_i)
+
+
+def cardinal(x, i, data_x):
+    return prod([(x - data_x[j]) / (data_x[i] - data_x[j]) for j in range(n) if i != j])
+
+
+def lagrange_polynomial(x, k, data_x, data_y):
+    return sum(data_y[i] * cardinal(x, i, data_x) for i in range(k))
+
+
+n = len(x_i)
+domain = (t, t_i[0], t_i[-1])
+
+# a)
+x_t = lagrange_polynomial(t, n, t_i, x_i)
+plot(
+    x_t,
+    domain,
+    title="x(t) – interpolated",
+    xlabel="t",
+    ylabel="x",
+    show=False,
+).save("x.png")
+
+# b)
+y_t = lagrange_polynomial(t, n, t_i, y_i)
+plot(
+    y_t,
+    domain,
+    title="y(t) – interpolated",
+    xlabel="t",
+    ylabel="y",
+    show=False,
+).save("y.png")
+
+# c) & d)
+interpolated = (x_t, y_t)
+actual = (cos(t) + sin(2 * t), sin(t))
+plot_parametric(
+    (x_t, y_t, domain, "interpolated"),
+    (cos(t) + sin(2 * t), sin(t), domain, "actual"),
+    title="helicopter path - approx. vs actual",
+    xlabel="x",
+    ylabel="y",
+    legend=True,
+    show=False,
+).save("path.png")