commit cdae52e238a65b7490caa996583df9fd4f3a9d58
parent 9419a6f6f31c19ef34666b1ac15371fc87fa9789
Author: miksa <milutin@popovic.xyz>
Date: Thu, 19 Jan 2023 12:52:48 +0000
adding nonlinear optimization
Diffstat:
13 files changed, 873 insertions(+), 0 deletions(-)
diff --git a/nlin_opt/build/sesh6.pdf b/nlin_opt/build/sesh6.pdf
Binary files differ.
diff --git a/nlin_opt/preamble.tex b/nlin_opt/preamble.tex
@@ -0,0 +1,58 @@
+\documentclass[a4paper]{article}
+
+\usepackage[T1]{fontenc}
+\usepackage[utf8]{inputenc}
+\usepackage{mlmodern}
+
+%\usepackage{ngerman} % Sprachanpassung Deutsch
+
+\usepackage{graphicx}
+\usepackage{geometry}
+\geometry{a4paper, top=15mm}
+
+\usepackage{subcaption}
+\usepackage[shortlabels]{enumitem}
+\usepackage{amssymb}
+\usepackage{amsthm}
+\usepackage{mathtools}
+\usepackage{braket}
+\usepackage{bbm}
+\usepackage{graphicx}
+\usepackage{float}
+\usepackage{yhmath}
+\usepackage{tikz}
+\usetikzlibrary{patterns,decorations.pathmorphing,positioning}
+\usetikzlibrary{calc,decorations.markings}
+
+%\usepackage[backend=biber, sorting=none]{biblatex}
+%\addbibresource{uni.bib}
+
+\usepackage[framemethod=TikZ]{mdframed}
+
+\tikzstyle{titlered} =
+ [draw=black, thick, fill=white,%
+ text=black, rectangle,
+ right, minimum height=.7cm]
+
+
+\usepackage[colorlinks=true,naturalnames=true,plainpages=false,pdfpagelabels=true]{hyperref}
+\usepackage[parfill]{parskip}
+\usepackage{lipsum}
+
+\usepackage{tcolorbox}
+\tcbuselibrary{skins,breakable}
+
+\pagestyle{myheadings}
+
+\newcommand{\eps}{\varepsilon}
+\usepackage[OT2,T1]{fontenc}
+\DeclareSymbolFont{cyrletters}{OT2}{wncyr}{m}{n}
+\DeclareMathSymbol{\Sha}{\mathalpha}{cyrletters}{"58}
+
+\markright{Popović\hfill Nonlinear Optimization\hfill}
+
+
+\title{University of Vienna\\ Faculty of Mathematics\\
+\vspace{1cm}Nonlinear Optimization Problems
+}
+\author{Milutin Popovic}
diff --git a/nlin_opt/prog/.ipynb_checkpoints/sesh6-checkpoint.ipynb b/nlin_opt/prog/.ipynb_checkpoints/sesh6-checkpoint.ipynb
@@ -0,0 +1,201 @@
+{
+ "cells": [
+ {
+ "cell_type": "code",
+ "execution_count": 1,
+ "id": "77f12d86",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "import numpy as np\n",
+ "import numdifftools as nd"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "921657fc",
+ "metadata": {},
+ "source": [
+ "Exercise 41:\n",
+ "-----------------\n",
+ "Implement the local Newton algorithm. Use as input data the starting vector x0,\n",
+ "the parameter for the stopping criterion ε and the parameter for the maximal\n",
+ "number of allowed iterations kmax. The sequence $x_0$, $x_1$, $x_2$, ...\n",
+ "containing the iteration history and the number of performed iterations should be returned.\n",
+ "The implemented algorithm should be tested for ε = 10^(−6) and kmax = 200, and the following\n",
+ "\n",
+ "Exercise 42:\n",
+ "-----------------\n",
+ "Implement the globalized Newton algorithm. Use as input data the starting vector $x_0$ the\n",
+ "parameter for the stopping criterion ε, the parameter for the maximal number of allowed\n",
+ "iterations kmax, the parameters for the determination of the Armijo step size σ and β, and\n",
+ "the parameters ρ and p. The sequence $x_0$, $x_1$, $x_2$, ... containing the iteration history and the\n",
+ "number of performed iterations should be returned.\n",
+ "The implemented algorithm should be tested for ε=10−6, kmax=200, ρ=10−8, p=2.1, σ=10−4"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 2,
+ "id": "bc388246",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "def local_newton(f, x_0, eps=10e-6, kmax=200):\n",
+ " f_grad = nd.Gradient(f)\n",
+ " f_hess = nd.Hessian(f)\n",
+ " x_k = x_0\n",
+ " \n",
+ " for k in range(kmax):\n",
+ " if np.linalg.norm(f_grad(x_k)) <= eps:\n",
+ " return x_k\n",
+ " break\n",
+ " \n",
+ " if type(x_0) == float:\n",
+ " d_k = -float(f_hess(x_k))/float(f_grad(x_k))\n",
+ " else:\n",
+ " d_k = np.linalg.solve(f_hess(x_k), -f_grad(x_k))\n",
+ " \n",
+ " x_k = x_k + np.array(d_k)\n",
+ " return x_k\n",
+ "\n",
+ "def global_newton(f, x_0, eps=10e-6, kmax=200, rho=10e-8, p=2.1, sig=10e-4, beta=0.5):\n",
+ " f_grad = nd.Gradient(f)\n",
+ " f_hess = nd.Hessian(f)\n",
+ " x_k = x_0\n",
+ " for k in range(kmax):\n",
+ " if np.linalg.norm(f_grad(x_k)) <= eps:\n",
+ " return x_k\n",
+ " break\n",
+ " \n",
+ " try:\n",
+ " if type(x_0) == float:\n",
+ " d_k = -float(f_hess(x_k))/float(f_grad(x_k))\n",
+ " else:\n",
+ " d_k = np.linalg.solve(f_hess(x_k), -f_grad(x_k))\n",
+ " if (np.dot(f_grad(x_k), d_k) > -rho*np.linalg.norm(d_k)**p):\n",
+ " d_k = -f_grad(x_k)\n",
+ " except np.linalg.LinAlgError:\n",
+ " d_k = -f_grad(x_k)\n",
+ " \n",
+ " # Find step size (Arnijno Rule)\n",
+ " t_k = 1\n",
+ " while f(x_k + t_k*d_k) > f(x_k) + sig*t_k * np.dot(f_grad(x_k), d_k):\n",
+ " t_k = t_k * beta\n",
+ " \n",
+ " x_k = x_k + t_k * d_k\n",
+ " \n",
+ " return x_k"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 3,
+ "id": "58828079",
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "a)\n",
+ "Local Newton Algorithm for x_0 = [-1.2, 1] is: [0.9999957 0.99999139]\n",
+ "Global Newton Algorithm for x_0 = [-1.2, 1] is: [1. 1.]\n",
+ "\n",
+ "b)\n",
+ "Local Newton Algorithm for x_0 = [0.25, 0.25, 0.25, 0.25] is: [1.50383723 0.7951105 0.47235406 0.32296666]\n",
+ "Global Newton Algorithm for x_0 = [0.25, 0.25, 0.25, 0.25] is: [1.50383723 0.7951105 0.47235406 0.32296666]\n",
+ "\n",
+ "c)\n",
+ "Local Newton Algorithm for x_0 = [1, 1] is: [158.75270666 79.36690168]\n",
+ "Global Newton Algorithm for x_0 = [1, 1] is: [1.24996948e-06 2.00000000e+06]\n",
+ "\n",
+ "d)\n",
+ "Local Newton Algorithm for x_0 = [-3, -1, -3, -1] is: [-1.41934291 2.02305705 0.36974654 0.12724277]\n",
+ "Global Newton Algorithm for x_0 = [-3, -1, -3, -1] is: [-1.41934305 2.02305745 0.36974706 0.12724316]\n",
+ "\n",
+ "e)\n",
+ "Local Newton Algorithm for x_0 = 2.0 is: -10.69705641570642\n",
+ "Local Newton Algorithm for x_0 = 1.0 is: 21.27903082235975\n",
+ "Local Newton Algorithm for x_0 = 0.5 is: 21.26959143622489\n",
+ "Global Newton Algorithm for x_0 = 2.0 is: -4.437802149414097e-11\n",
+ "Global Newton Algorithm for x_0 = 1.0 is: 8.446339466597688e-14\n",
+ "Global Newton Algorithm for x_0 = 0.5 is: 7.344043315106116e-14\n",
+ "\n"
+ ]
+ }
+ ],
+ "source": [
+ "f_a = lambda x: (1-x[0])**2 + 100*(x[1] - x[0]**2)**2\n",
+ "x_a_0 = [-1.2, 1]\n",
+ "print('a)')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_a_0} is: {local_newton(f_a, x_a_0)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_a_0} is: {global_newton(f_a, x_a_0)}')\n",
+ "print('')\n",
+ "\n",
+ "f_b = lambda x: sum( (sum(np.cos(x[j]) + (i+1)*(1 - np.cos(x[i])) - np.sin(x[i]) for j in range(4) ) )**2 for i in range(4) )\n",
+ "x_b_0 = [1/4, 1/4, 1/4, 1/4]\n",
+ "print('b)')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_b_0} is: {local_newton(f_b, x_b_0)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_b_0} is: {global_newton(f_b, x_b_0)}')\n",
+ "print('')\n",
+ "\n",
+ "f_c = lambda x: (x[0] - 10**6)**2 + (x[1] - 2*10**6)**2 + (x[0]*x[1] - 2)**2\n",
+ "x_c_0 = [1, 1]\n",
+ "print('c)')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_c_0} is: {local_newton(f_c, x_c_0)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_c_0} is: {global_newton(f_c, x_c_0)}')\n",
+ "print('')\n",
+ "\n",
+ "f_d = lambda x: 100*(x[1] - x[0]**2)**2 + (1-x[0])**2 + 90*(x[3] - x[2]**2)**2 + (1 - x[2])**2 + 10*(x[1] - x[3] - 2)**2 + 1/10*(x[1] - x[3])**2\n",
+ "x_d_0 = [-3, -1, -3, -1]\n",
+ "print('d)')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_d_0} is: {local_newton(f_d, x_d_0)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_d_0} is: {global_newton(f_d, x_d_0)}')\n",
+ "print('')\n",
+ "\n",
+ "f_e = lambda x: np.sqrt(1 + x**2)\n",
+ "x_e_00 = 2.0\n",
+ "x_e_01 = 1.0\n",
+ "x_e_02 = 1/2\n",
+ "print('e)')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_e_00} is: {local_newton(f_e, x_e_00)}')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_e_01} is: {local_newton(f_e, x_e_01)}')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_e_02} is: {local_newton(f_e, x_e_02)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_e_00} is: {global_newton(f_e, x_e_00)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_e_01} is: {global_newton(f_e, x_e_01)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_e_02} is: {global_newton(f_e, x_e_02)}')\n",
+ "print('')"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "b1f45770",
+ "metadata": {},
+ "outputs": [],
+ "source": []
+ }
+ ],
+ "metadata": {
+ "kernelspec": {
+ "display_name": "Python 3",
+ "language": "python",
+ "name": "python3"
+ },
+ "language_info": {
+ "codemirror_mode": {
+ "name": "ipython",
+ "version": 3
+ },
+ "file_extension": ".py",
+ "mimetype": "text/x-python",
+ "name": "python",
+ "nbconvert_exporter": "python",
+ "pygments_lexer": "ipython3",
+ "version": "3.10.9"
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}
diff --git a/nlin_opt/prog/.ipynb_checkpoints/shesh6-checkpoint.ipynb b/nlin_opt/prog/.ipynb_checkpoints/shesh6-checkpoint.ipynb
@@ -0,0 +1,192 @@
+{
+ "cells": [
+ {
+ "cell_type": "code",
+ "execution_count": 66,
+ "id": "d441cf6d",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "import numpy as np\n",
+ "import numdifftools as nd"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "b2a1b9b1",
+ "metadata": {},
+ "source": [
+ "Exercise 41:\n",
+ "-----------------\n",
+ "Implement the local Newton algorithm. Use as input data the starting vector x0,\n",
+ "the parameter for the stopping criterion ε and the parameter for the maximal\n",
+ "number of allowed iterations kmax. The sequence $x_0$, $x_1$, $x_2$, ...\n",
+ "containing the iteration history and the number of performed iterations should be returned.\n",
+ "The implemented algorithm should be tested for ε = 10^(−6) and kmax = 200, and the following\n",
+ "\n",
+ "Exercise 42:\n",
+ "-----------------\n",
+ "Implement the globalized Newton algorithm. Use as input data the starting vector $x_0$ the\n",
+ "parameter for the stopping criterion ε, the parameter for the maximal number of allowed\n",
+ "iterations kmax, the parameters for the determination of the Armijo step size σ and β, and\n",
+ "the parameters ρ and p. The sequence $x_0$, $x_1$, $x_2$, ... containing the iteration history and the\n",
+ "number of performed iterations should be returned.\n",
+ "The implemented algorithm should be tested for ε=10−6, kmax=200, ρ=10−8, p=2.1, σ=10−4"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 257,
+ "id": "89869b1f",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "def local_newton(f, x_0, eps=10e-6, kmax=200):\n",
+ " f_grad = nd.Gradient(f)\n",
+ " f_hess = nd.Hessian(f)\n",
+ " x_k = x_0\n",
+ " \n",
+ " for k in range(kmax):\n",
+ " if np.linalg.norm(f_grad(x_k)) <= eps:\n",
+ " return x_k\n",
+ " break\n",
+ " \n",
+ " if type(x_0) == float:\n",
+ " d_k = -float(f_hess(x_k))/float(f_grad(x_k))\n",
+ " else:\n",
+ " d_k = np.linalg.solve(f_hess(x_k), -f_grad(x_k))\n",
+ " \n",
+ " x_k = x_k + np.array(d_k)\n",
+ " return x_k\n",
+ "\n",
+ "def global_newton(f, x_0, eps=10e-6, kmax=200, rho=10e-8, p=2.1, sig=10e-4, beta=0.5):\n",
+ " f_grad = nd.Gradient(f)\n",
+ " f_hess = nd.Hessian(f)\n",
+ " x_k = x_0\n",
+ " for k in range(kmax):\n",
+ " if np.linalg.norm(f_grad(x_k)) <= eps:\n",
+ " return x_k\n",
+ " break\n",
+ " \n",
+ " try:\n",
+ " if type(x_0) == float:\n",
+ " d_k = -float(f_hess(x_k))/float(f_grad(x_k))\n",
+ " else:\n",
+ " d_k = np.linalg.solve(f_hess(x_k), -f_grad(x_k))\n",
+ " # if the linear system above has no solutino or a condition is satisfied\n",
+ " except np.linalg.LinAlgError or (np.dot(f_grad(x_k), d_k) > -rho*np.linalg.norm(d_k)**p):\n",
+ " d_k = -f_grad(x_k)\n",
+ " \n",
+ " # Find step size (Arnijno Rule)\n",
+ " t_k = 1\n",
+ " while f(x_k + t_k*d_k) > f(x_k) + sig*t_k * np.dot(f_grad(x_k), d_k):\n",
+ " t_k = t_k * beta\n",
+ " \n",
+ " x_k = x_k + t_k * d_k\n",
+ " \n",
+ " return x_k"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 266,
+ "id": "0f5d7d60",
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "a)\n",
+ "Local Newton Algorithm for x_0 = [-1.2, 1] is: [0.9999957 0.99999139]\n",
+ "Global Newton Algorithm for x_0 = [-1.2, 1] is: [1. 1.]\n",
+ "\n",
+ "b)\n",
+ "Local Newton Algorithm for x_0 = [0.25, 0.25, 0.25, 0.25] is: [1.50383723 0.7951105 0.47235406 0.32296666]\n",
+ "Global Newton Algorithm for x_0 = [0.25, 0.25, 0.25, 0.25] is: [1.50383723 0.7951105 0.47235406 0.32296666]\n",
+ "\n",
+ "c)\n",
+ "Local Newton Algorithm for x_0 = [1, 1] is: [158.75270666 79.36690168]\n",
+ "Global Newton Algorithm for x_0 = [1, 1] is: [120.72507856 141.61646967]\n",
+ "\n",
+ "d)\n",
+ "Local Newton Algorithm for x_0 = [-3, -1, -3, -1] is: [-1.41934291 2.02305705 0.36974654 0.12724277]\n",
+ "Global Newton Algorithm for x_0 = [-3, -1, -3, -1] is: [-1.41934305 2.02305745 0.36974706 0.12724316]\n",
+ "\n",
+ "e)\n",
+ "Local Newton Algorithm for x_0 = 2.0 is: -10.69705641570642\n",
+ "Local Newton Algorithm for x_0 = 1.0 is: 21.27903082235975\n",
+ "Local Newton Algorithm for x_0 = 0.5 is: 21.26959143622489\n",
+ "Global Newton Algorithm for x_0 = 2.0 is: -5.065135369566478e-06\n",
+ "Global Newton Algorithm for x_0 = 1.0 is: -7.786539613758756e-06\n",
+ "Global Newton Algorithm for x_0 = 0.5 is: -7.832383476097789e-06\n",
+ "\n"
+ ]
+ }
+ ],
+ "source": [
+ "f_a = lambda x: (1-x[0])**2 + 100*(x[1] - x[0]**2)**2\n",
+ "x_a_0 = [-1.2, 1]\n",
+ "print('a)')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_a_0} is: {local_newton(f_a, x_a_0)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_a_0} is: {global_newton(f_a, x_a_0)}')\n",
+ "print('')\n",
+ "\n",
+ "f_b = lambda x: sum( (sum(np.cos(x[j]) + (i+1)*(1 - np.cos(x[i])) - np.sin(x[i]) for j in range(4) ) )**2 for i in range(4) )\n",
+ "x_b_0 = [1/4, 1/4, 1/4, 1/4]\n",
+ "print('b)')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_b_0} is: {local_newton(f_b, x_b_0)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_b_0} is: {global_newton(f_b, x_b_0)}')\n",
+ "print('')\n",
+ "\n",
+ "f_c = lambda x: (x[0] - 10**6)**2 + (x[1] - 2*10**6)**2 + (x[0]*x[1] - 2)**2\n",
+ "x_c_0 = [1, 1]\n",
+ "print('c)')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_c_0} is: {local_newton(f_c, x_c_0)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_c_0} is: {global_newton(f_c, x_c_0)}')\n",
+ "print('')\n",
+ "\n",
+ "f_d = lambda x: 100*(x[1] - x[0]**2)**2 + (1-x[0])**2 + 90*(x[3] - x[2]**2)**2 + (1 - x[2])**2 + 10*(x[1] - x[3] - 2)**2 + 1/10*(x[1] - x[3])**2\n",
+ "x_d_0 = [-3, -1, -3, -1]\n",
+ "print('d)')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_d_0} is: {local_newton(f_d, x_d_0)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_d_0} is: {global_newton(f_d, x_d_0)}')\n",
+ "print('')\n",
+ "\n",
+ "f_e = lambda x: np.sqrt(1 + x**2)\n",
+ "x_e_00 = 2.0\n",
+ "x_e_01 = 1.0\n",
+ "x_e_02 = 1/2\n",
+ "print('e)')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_e_00} is: {local_newton(f_e, x_e_00)}')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_e_01} is: {local_newton(f_e, x_e_01)}')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_e_02} is: {local_newton(f_e, x_e_02)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_e_00} is: {global_newton(f_e, x_e_00)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_e_01} is: {global_newton(f_e, x_e_01)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_e_02} is: {global_newton(f_e, x_e_02)}')\n",
+ "print('')"
+ ]
+ }
+ ],
+ "metadata": {
+ "kernelspec": {
+ "display_name": "Python 3",
+ "language": "python",
+ "name": "python3"
+ },
+ "language_info": {
+ "codemirror_mode": {
+ "name": "ipython",
+ "version": 3
+ },
+ "file_extension": ".py",
+ "mimetype": "text/x-python",
+ "name": "python",
+ "nbconvert_exporter": "python",
+ "pygments_lexer": "ipython3",
+ "version": "3.10.9"
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}
diff --git a/nlin_opt/prog/sesh6.ipynb b/nlin_opt/prog/sesh6.ipynb
@@ -0,0 +1,209 @@
+{
+ "cells": [
+ {
+ "cell_type": "markdown",
+ "id": "bcc85c40",
+ "metadata": {},
+ "source": [
+ "# MILUTIN POPOVIC: EXERCISE SHEET 6"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 1,
+ "id": "77f12d86",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "import numpy as np\n",
+ "import numdifftools as nd"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "921657fc",
+ "metadata": {},
+ "source": [
+ "Exercise 41:\n",
+ "-----------------\n",
+ "Implement the local Newton algorithm. Use as input data the starting vector x0,\n",
+ "the parameter for the stopping criterion ε and the parameter for the maximal\n",
+ "number of allowed iterations kmax. The sequence $x_0$, $x_1$, $x_2$, ...\n",
+ "containing the iteration history and the number of performed iterations should be returned.\n",
+ "The implemented algorithm should be tested for ε = 10^(−6) and kmax = 200, and the following\n",
+ "\n",
+ "Exercise 42:\n",
+ "-----------------\n",
+ "Implement the globalized Newton algorithm. Use as input data the starting vector $x_0$ the\n",
+ "parameter for the stopping criterion ε, the parameter for the maximal number of allowed\n",
+ "iterations kmax, the parameters for the determination of the Armijo step size σ and β, and\n",
+ "the parameters ρ and p. The sequence $x_0$, $x_1$, $x_2$, ... containing the iteration history and the\n",
+ "number of performed iterations should be returned.\n",
+ "The implemented algorithm should be tested for ε=10−6, kmax=200, ρ=10−8, p=2.1, σ=10−4"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 2,
+ "id": "bc388246",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "def local_newton(f, x_0, eps=10e-6, kmax=200):\n",
+ " f_grad = nd.Gradient(f)\n",
+ " f_hess = nd.Hessian(f)\n",
+ " x_k = x_0\n",
+ " \n",
+ " for k in range(kmax):\n",
+ " if np.linalg.norm(f_grad(x_k)) <= eps:\n",
+ " return x_k\n",
+ " break\n",
+ " \n",
+ " if type(x_0) == float:\n",
+ " d_k = -float(f_hess(x_k))/float(f_grad(x_k))\n",
+ " else:\n",
+ " d_k = np.linalg.solve(f_hess(x_k), -f_grad(x_k))\n",
+ " \n",
+ " x_k = x_k + np.array(d_k)\n",
+ " return x_k\n",
+ "\n",
+ "def global_newton(f, x_0, eps=10e-6, kmax=200, rho=10e-8, p=2.1, sig=10e-4, beta=0.5):\n",
+ " f_grad = nd.Gradient(f)\n",
+ " f_hess = nd.Hessian(f)\n",
+ " x_k = x_0\n",
+ " for k in range(kmax):\n",
+ " if np.linalg.norm(f_grad(x_k)) <= eps:\n",
+ " return x_k\n",
+ " break\n",
+ " \n",
+ " try:\n",
+ " if type(x_0) == float:\n",
+ " d_k = -float(f_hess(x_k))/float(f_grad(x_k))\n",
+ " else:\n",
+ " d_k = np.linalg.solve(f_hess(x_k), -f_grad(x_k))\n",
+ " if (np.dot(f_grad(x_k), d_k) > -rho*np.linalg.norm(d_k)**p):\n",
+ " d_k = -f_grad(x_k)\n",
+ " except np.linalg.LinAlgError:\n",
+ " d_k = -f_grad(x_k)\n",
+ " \n",
+ " # Find step size (Arnijno Rule)\n",
+ " t_k = 1\n",
+ " while f(x_k + t_k*d_k) > f(x_k) + sig*t_k * np.dot(f_grad(x_k), d_k):\n",
+ " t_k = t_k * beta\n",
+ " \n",
+ " x_k = x_k + t_k * d_k\n",
+ " \n",
+ " return x_k"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 3,
+ "id": "58828079",
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "a)\n",
+ "Local Newton Algorithm for x_0 = [-1.2, 1] is: [0.9999957 0.99999139]\n",
+ "Global Newton Algorithm for x_0 = [-1.2, 1] is: [1. 1.]\n",
+ "\n",
+ "b)\n",
+ "Local Newton Algorithm for x_0 = [0.25, 0.25, 0.25, 0.25] is: [1.50383723 0.7951105 0.47235406 0.32296666]\n",
+ "Global Newton Algorithm for x_0 = [0.25, 0.25, 0.25, 0.25] is: [1.50383723 0.7951105 0.47235406 0.32296666]\n",
+ "\n",
+ "c)\n",
+ "Local Newton Algorithm for x_0 = [1, 1] is: [158.75270666 79.36690168]\n",
+ "Global Newton Algorithm for x_0 = [1, 1] is: [1.24996948e-06 2.00000000e+06]\n",
+ "\n",
+ "d)\n",
+ "Local Newton Algorithm for x_0 = [-3, -1, -3, -1] is: [-1.41934291 2.02305705 0.36974654 0.12724277]\n",
+ "Global Newton Algorithm for x_0 = [-3, -1, -3, -1] is: [-1.41934305 2.02305745 0.36974706 0.12724316]\n",
+ "\n",
+ "e)\n",
+ "Local Newton Algorithm for x_0 = 2.0 is: -10.69705641570642\n",
+ "Local Newton Algorithm for x_0 = 1.0 is: 21.27903082235975\n",
+ "Local Newton Algorithm for x_0 = 0.5 is: 21.26959143622489\n",
+ "Global Newton Algorithm for x_0 = 2.0 is: -4.437802149414097e-11\n",
+ "Global Newton Algorithm for x_0 = 1.0 is: 8.446339466597688e-14\n",
+ "Global Newton Algorithm for x_0 = 0.5 is: 7.344043315106116e-14\n",
+ "\n"
+ ]
+ }
+ ],
+ "source": [
+ "f_a = lambda x: (1-x[0])**2 + 100*(x[1] - x[0]**2)**2\n",
+ "x_a_0 = [-1.2, 1]\n",
+ "print('a)')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_a_0} is: {local_newton(f_a, x_a_0)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_a_0} is: {global_newton(f_a, x_a_0)}')\n",
+ "print('')\n",
+ "\n",
+ "f_b = lambda x: sum( (sum(np.cos(x[j]) + (i+1)*(1 - np.cos(x[i])) - np.sin(x[i]) for j in range(4) ) )**2 for i in range(4) )\n",
+ "x_b_0 = [1/4, 1/4, 1/4, 1/4]\n",
+ "print('b)')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_b_0} is: {local_newton(f_b, x_b_0)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_b_0} is: {global_newton(f_b, x_b_0)}')\n",
+ "print('')\n",
+ "\n",
+ "f_c = lambda x: (x[0] - 10**6)**2 + (x[1] - 2*10**6)**2 + (x[0]*x[1] - 2)**2\n",
+ "x_c_0 = [1, 1]\n",
+ "print('c)')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_c_0} is: {local_newton(f_c, x_c_0)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_c_0} is: {global_newton(f_c, x_c_0)}')\n",
+ "print('')\n",
+ "\n",
+ "f_d = lambda x: 100*(x[1] - x[0]**2)**2 + (1-x[0])**2 + 90*(x[3] - x[2]**2)**2 + (1 - x[2])**2 + 10*(x[1] - x[3] - 2)**2 + 1/10*(x[1] - x[3])**2\n",
+ "x_d_0 = [-3, -1, -3, -1]\n",
+ "print('d)')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_d_0} is: {local_newton(f_d, x_d_0)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_d_0} is: {global_newton(f_d, x_d_0)}')\n",
+ "print('')\n",
+ "\n",
+ "f_e = lambda x: np.sqrt(1 + x**2)\n",
+ "x_e_00 = 2.0\n",
+ "x_e_01 = 1.0\n",
+ "x_e_02 = 1/2\n",
+ "print('e)')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_e_00} is: {local_newton(f_e, x_e_00)}')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_e_01} is: {local_newton(f_e, x_e_01)}')\n",
+ "print(f'Local Newton Algorithm for x_0 = {x_e_02} is: {local_newton(f_e, x_e_02)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_e_00} is: {global_newton(f_e, x_e_00)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_e_01} is: {global_newton(f_e, x_e_01)}')\n",
+ "print(f'Global Newton Algorithm for x_0 = {x_e_02} is: {global_newton(f_e, x_e_02)}')\n",
+ "print('')"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "b1f45770",
+ "metadata": {},
+ "outputs": [],
+ "source": []
+ }
+ ],
+ "metadata": {
+ "kernelspec": {
+ "display_name": "Python 3",
+ "language": "python",
+ "name": "python3"
+ },
+ "language_info": {
+ "codemirror_mode": {
+ "name": "ipython",
+ "version": 3
+ },
+ "file_extension": ".py",
+ "mimetype": "text/x-python",
+ "name": "python",
+ "nbconvert_exporter": "python",
+ "pygments_lexer": "ipython3",
+ "version": "3.10.9"
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}
diff --git a/nlin_opt/sesh6.tex b/nlin_opt/sesh6.tex
@@ -0,0 +1,213 @@
+\include{./preamble.tex}
+
+\begin{document}
+\maketitle
+\tableofcontents
+\section{Sheet 6}
+\subsection{Exercise 35}
+Let $M \in \mathbb{R}^{n \times n}$ with $\|M\| < 1$. Show that $I - M$ is
+regular and
+\begin{equation}
+ \|\left( I - M \right) \| \le \frac{1}{1-\|M\|}. \label{eq: ineq}
+\end{equation}
+Suppose $I-M$ is not singular then for $x \in \mathbb{R}^{n}$ we have that
+\begin{align}
+ &\left( I - M \right)x = 0\\
+ \Leftrightarrow &Ix - Mx = 0\\
+ \Leftrightarrow & Mx = x.
+\end{align}
+But since $\|M\|<1$ then $\forall x \in \mathbb{R}^{n}$ we have that $\|Mx\|
+< \|x\|$. This means that
+\begin{align}
+ \text{ker}\left(I-M \right) = \emptyset,
+\end{align}
+so $I-M$ is regular. The identity on the other hand is derived by the
+following observation
+\begin{align}
+ \left( I-M \right)^{-1} - \left( I-M \right)^{-1} M = \left( I -M
+ \right)\left( I-M \right)^{-1} = I,
+\end{align}
+Then we calculate
+\begin{align}
+ \|\left(I-M \right)^{-1}\|
+ &= \|I+\left( I-M \right)^{-1} M\|\\
+ &\le \|I\| + \|(I-M)^{-1}\|\|M\|,
+\end{align}
+rearranging gives
+\begin{align}
+ \|\left( I-M \right)^{-1}\| - \|(I-M)^{-1}\| \le \|I\|\\
+ \|(I-M)^{-1}\|\left(1-\|M\| \right) \le 1\\
+ \|\left( I-M \right)^{-1}\| \le \frac{1}{1-\|M\|}.
+\end{align}
+Now let $A,B \in \mathbb{R}^{n \times n}$ with $\|I-BA\|<1$. Show that $A$
+and $B$ are regular and that
+\begin{align}
+ \|B^{-1}\| \le \frac{\|A\|}{1-\|I-BA\|}\\
+ \|A^{-1}\| \le \frac{\|B\|}{1-\|I-BA\|}\\
+\end{align}
+We know that for $M \in \mathbb{R}^{n \times n}$ with $\|M\|<1$ then $I-M$
+is regular and the inequality in \ref{eq: ineq} holds. Set $M =
+I-BA$ then $I-M=AB$ is regular. Because $AB$ is regular so is $A$ and $B$.
+Now note that for all regular matrices we have that $\|A^{-1}\|\le
+\|A\|^{-1}$. Furthermore
+\begin{align}
+ \|\left(AB \right)^{-1}\| \le \|B^{-1}\|\|A^{-1}\|.
+\end{align}
+Then for $A$ we have
+\begin{align}
+ \|A^{-1}\|\le \frac{1}{\|B^{-1}\|}\frac{1}{1-\|I-BA\|}\le
+ \frac{\|B\|}{1-\|I-BA\|}.
+\end{align}
+and for $B$
+\begin{align}
+ \|B^{-1}\|\le \frac{1}{\|A^{-1}\|}\frac{1}{1-\|I-BA\|}\le
+ \frac{\|A\|}{1-\|I-BA\|}.
+\end{align}
+\subsection{Exercise 36}
+Let $f:\mathbb{R}^{2}\to \mathbb{R}$, $f(x, y) = x^4 + 2x^2y^2 + y^4$. Show
+that the local Newton algorithm converges to the unique global minimum of
+$f$ for every $(x^0, y^0) \in \mathbb{R}^{2}\backslash\{(0, 0)^T\}$.
+First we determine the minimum $x^*$. Note that $f(x, y) = \left( x^2 + y^2
+\right)^2 \ge 0$ for all $(x, y)^T \in \mathbb{R}^{2}$. Since $f$ is strongly convex the
+only minimum, which is the global minimum is $(x, y)^T = (0 ,0)^T$. The
+Hessian of $f$ is
+\begin{align}
+ \nabla^2 f(x,y) = \begin{pmatrix}12x^2+4y^2 & 8xy\\ 8xy & 12y^2 + 4x^2
+ \end{pmatrix}.
+\end{align}
+Now note that the Hessian at the minimum $\nabla^2f(0,0)$ is the zero matrix
+which is singular. But considering starting vectors $(x, y)^T \neq (0,
+0)^T$, all we need in the local Newton algorithm is the solution to the equation
+$\nabla^2f(x^k)d_k = -\nabla f(x^k)$. Meaning we need to show that
+$\nabla^2f(x,y)$ is regular for all $(x, y)^T \neq (0,0)^T$, in this
+case we look at the determinant of the Hessian
+\begin{align}
+ \det(\nabla^2 f(x,y)) = 48*(x^2 + y^2)^2 > 0 \;\; \forall (x, y)^T \neq (0,
+ 0)^T.
+\end{align}
+This means that the sequence $(x^k)_{k\ge 0}$ produced by the local newton
+algorithm converges to the unique global minimum of $f$ given by $(0, 0)^T$.
+Indeed if we calculate the solution of the system $\nabla^2f(x,y)d =
+-\nabla f(x,y)$ we get that $d = (\frac{x}{3}, \frac{y}{3})^T$.
+
+\subsection{Exercise 37}
+Show that the local Newton Algorithm is invariant to affine-linear
+transformation, for a regular matrix $A \in \mathbb{R}^{n \times n}$ and $c
+\in \mathbb{R}^{n}$, $(x^k)_{k\ge 0}$ the sequence generated by the local
+Newton algorithm for minimizing $f$ with starting vector $x^0$. Then let
+$(y^k)_{k\ge 0}$ the sequence generated by the local Newton algorithm for the
+function $g(y) := f(Ay +c)$ with starting vector $y^0$, then
+\begin{align}
+ x^0 = Ay^0 + c \;\; \implies x^k = Ay^k + c \;\;\; \forall k \ge 0.
+\end{align}
+First of all we calculate the gradient and the hessian for $g$
+\begin{align}
+ \nabla g(y) = \nabla f(Ay+c) = A^T \nabla f(Ay+c)\\
+ \nabla^2 g(y) = \nabla^2 f(Ay+c) = A^T \nabla^2 f(Ay+c) A\\
+\end{align}
+Now we need to prove that $x^{k+1} = Ay^{k+1} + c$
+\begin{align}
+ x^{k+1}q
+ &= x^{k} + d_k = x^{k} - \left( \nabla^2 f(x^{k})\right)^{-1}
+ \nabla f(x_k)\\
+ &=Ay^k + c - \left( \nabla^2(f(Ay^{k}+c )\right)^{-1} \nabla f(Ay^{k}+c)\\
+ &=Ay^{k} + c - A A^{-1}\left( \nabla^2(f(Ay^{k}+c )\right)^{-1} \nabla f(Ay^{k}+c)\\
+ &=Ay^{k} + c - A A^{-1} A^T A^{-T}\left( \nabla^2(f(Ay^{k}+c )\right)^{-1} \nabla f(Ay^{k}+c)\\
+ &=A\left(y^{k} - \left( A^{T}\nabla^2 f(Ay^{k}+c)A\right)^{-1}A^T \nabla
+ f(Ay^{k}+c ) \right) +c \\
+ &= Ay^{k+1} +c.
+\end{align}
+by that the induction is finished.
+\subsection{Exercise 38}
+Let $M \in \mathbb{R}^{n \times n}$ be a regular matrix and $\{M\}_{k\ge 0}
+\in \mathbb{R}^{n \times n}$ a sequence of matrices which converge to M as $k
+\to \infty $. Sow that there exists a $k_0 \ge 0$ such that $M_k$ is regular
+for all $k\ge k_0$ and the sequence $\{M_k^{-1}\}_{k\ge 0}$ converges to
+$M^{-1}$.
+\newline
+The map $M \to M^{-1}$ is a continuous invertible meaning it is monotone.
+$M^{-1} = \frac{\text{adj}(M)}{\text{det}(M)}$. Then convergence means that
+there is a $k\ge k_0$ such that for all $M_k \in B_{\frac{1}{k}}(M)$ we have
+that $\|M_k -M\| < \frac{1}{k}$ then $M_k$ is sufficiently close to $M$ and
+so regular. Since $\{M_k\}_{k\ge k_0} \cup {M}$ is a compact set
+of invertible matrices so is $\{M_k^{-1}\}_{k\ge k_0} \cup {M^{-1}}$,
+meaning it is bounded. This means that $\{M^{-1}_k\}_{k\ge k_0}$ converges to
+$M^{-1}$.
+\subsection{Exercise 39}
+Let $H \in \mathbb{R}^{n \times n}$ be regular $u, v \in \mathbb{R}^{n}$
+arbitrary. Show that $H+uv^T$ regular $\Leftrightarrow$ $1+v^TH^{-1}u \neq
+0$, then the Sherman-Morrison formula holds
+\begin{align}
+ \left( H + uv^T \right)^{-1} = \left( I - \frac{1}{1-v^TH^{-1} u} H^{-1}
+ uv^T\right)H^{-1} \label{eq: smf}
+\end{align}
+Let $1+ v^TH^{-1}u =0$ then
+\begin{align}
+ \text{det}\left( H + uv^T \right) =(1+v^T H^{-1} u)\text{det}(H) = 0.
+\end{align}
+This means that $H$ is not invertible. Now we need to check if the inverse
+really holds which is done by simply multiplying
+\begin{align}
+ &\left( H+uv^{T} \right) \left( H^{-1} - \frac{H^{-1} uv^T
+ H^{-1}}{1+v^{T}H^{-1}u} \right) = \\
+ &=H H^{-1} + uv^T H^{-1} - H \frac{H^{-1}uv^T H^{-1}}{1+v^TH^{-1}u}
+ uv^{T}\frac{H^{-1}uv^TH^{-1}}{1+v^{T}H^{-1}u}\\
+ &=I + uv^T H^{-1} - \frac{uv^T H^{-1} + uv^T H^{-1}uv^T
+ H^{-1}}{1+v^TH^{-1}u}\\
+ &=I+uv^TH^{-1} - \frac{u\left( 1+v^{T}H^{-1}u \right)
+ v^{T}H^{-1}}{1+v^{T}H^{-1} u}\\
+ &=I + uv^{T}H^{-1} - uv^{T}H^{-1} \\
+ &= I
+\end{align}
+Since these are square matrices $AB=I$ is the same as $BA=I$.
+\subsection{Exercise 40}
+Consider the quadratic optimization problem
+\begin{align}
+ \min \quad &f(x):=\gamma + c^{T}x + \frac{1}{2}x^{T}Qx, \label{eq: opp}\\
+ \text{s.t}\quad &h(x) := b^{T}x = 0, \nonumber
+\end{align}
+with $Q \in \mathbb{R}^{n \times n}$ SPD , $b, c \in \mathbb{R}^{n}$, $b\neq
+0$ and $\gamma \in \mathbb{R}$. For a given $\alpha > 0$ find the minimum
+$x^*(\alpha)$ of the penalty function
+\begin{align}
+ P(x;\alpha) := f(x) + \frac{\alpha}{2}\left( h(x) \right)^{2}
+\end{align}
+determine $x^* := \lim_{\alpha \to \infty}x^{*}(\alpha)$ and prove that
+$x^{*}$ is a unique optimal solution of the optimization problem in \ref{eq:
+opp}. We start with calculating the minimum of $P(x(\alpha))$
+\begin{align}
+ \nabla P(x(\alpha))
+ &= \nabla f(x) + \frac{\alpha}{2}\nabla h(x)^{2}\\
+ &= c + Qx + \frac{\alpha}{2} 2 h(x) \nabla h(x) \\
+ &=c + Qx + \alpha b^T x b = 0\\
+ \nonumber\\
+ Qx + \alpha b b ^T x = -c \\
+ \left( Q + \alpha b b^T \right) x = -c.
+\end{align}
+Using the Sherman-Morrison formula in \ref{eq: smf} for $H = Q$, $u=\alpha
+b$ and $v = b$ we get
+\begin{align}
+ x^{*}(\alpha) = \left( \frac{\alpha}{1+\alpha b^{T}Q^{-1}b}Q^{-1} bb^{T}
+ -I\right) Q^{-1} c
+\end{align}
+The limit is then (the standard limit $\frac{x}{1+kx} \to \frac{1}{k}$ as $x$
+goes to infinity)
+\begin{align}
+ x^{*} &= \lim_{\alpha \to \infty}x^{*}(\alpha)\\
+ &= \left(\frac{Q^{-1}bb^{T}}{b^{T}Q^{-1}b} -I \right) Q^{-1} c.
+\end{align}
+To show that $x^{*}$ is a unique solution of the optimization problem
+\ref{eq: opp} we need to show that it satisfies the optimality condition and
+that it is unique. Now it is unique because $Q$ is SPD meaning it is regular
+and invertible and $\nabla^2 f = Q>0$. Further more $(x^{*}, \alpha)$ is a
+KKT point of $P(x, \alpha)$ then $x^{*}$ is a minumum of the optimization
+problem. Now we show that $b^{T}x^{*} = 0$:
+\begin{align}
+ b^{T}x^{*} &= \left(\frac{b^{T}Q^{-1}bb^{T}Q^{-1}}{b^{T}Q^{-1}b}-
+ b^{T}Q^{-1} \right) c \\
+ &= \left( b^{T}Q - b^{T}Q \right) c \\
+ &= 0.
+\end{align}
+\end{document}
+
+
diff --git a/nlin_opt/sheets/Nonlinear_optimization skript.pdf b/nlin_opt/sheets/Nonlinear_optimization skript.pdf
Binary files differ.
diff --git a/nlin_opt/sheets/exercisesession1.pdf b/nlin_opt/sheets/exercisesession1.pdf
Binary files differ.
diff --git a/nlin_opt/sheets/exercisesession2.pdf b/nlin_opt/sheets/exercisesession2.pdf
Binary files differ.
diff --git a/nlin_opt/sheets/exercisesession3.pdf b/nlin_opt/sheets/exercisesession3.pdf
Binary files differ.
diff --git a/nlin_opt/sheets/exercisesession4.pdf b/nlin_opt/sheets/exercisesession4.pdf
Binary files differ.
diff --git a/nlin_opt/sheets/exercisesession5.pdf b/nlin_opt/sheets/exercisesession5.pdf
Binary files differ.
diff --git a/nlin_opt/sheets/exercisesession6.pdf b/nlin_opt/sheets/exercisesession6.pdf
Binary files differ.
