diff --git a/README.md b/README.md index 8059993..73ca14a 100644 --- a/README.md +++ b/README.md @@ -74,8 +74,8 @@ general, I will not post homework solutions. |3|1/31||Consistent Coding habits| | |2/2|5|Root Finding| |4|2/7|6|Root Finding con’d| -| |2/9|7|Optimization| -|5|2/14||Intro to Linear Algebra| +| |2/9|7| **Snow Day**| +|5|2/14|| Optimization | | |2/16|8|Linear Algebra| |6|2/21|9|Linear systems: Gauss elimination| | |2/23|10|Linear Systems: LU factorization| diff --git a/lecture_07/.lecture_07.md.swo b/lecture_07/.lecture_07.md.swo deleted file mode 100644 index 244c92b..0000000 Binary files a/lecture_07/.lecture_07.md.swo and /dev/null differ diff --git a/lecture_08/Auchain_model.png b/lecture_08/Auchain_model.png new file mode 100644 index 0000000..dce9cd7 Binary files /dev/null and b/lecture_08/Auchain_model.png differ diff --git a/lecture_08/goldmin.m b/lecture_08/goldmin.m new file mode 100644 index 0000000..fb4ce0b --- /dev/null +++ b/lecture_08/goldmin.m @@ -0,0 +1,36 @@ +function [x,fx,ea,iter]=goldmin(f,xl,xu,es,maxit,varargin) +% goldmin: minimization golden section search +% [x,fx,ea,iter]=goldmin(f,xl,xu,es,maxit,p1,p2,...): +% uses golden section search to find the minimum of f +% input: +% f = name of function +% xl, xu = lower and upper guesses +% es = desired relative error (default = 0.0001%) +% maxit = maximum allowable iterations (default = 50) +% p1,p2,... = additional parameters used by f +% output: +% x = location of minimum +% fx = minimum function value +% ea = approximate relative error (%) +% iter = number of iterations +if nargin<3,error('at least 3 input arguments required'),end +if nargin<4|isempty(es), es=0.0001;end +if nargin<5|isempty(maxit), maxit=50;end +phi=(1+sqrt(5))/2; +iter=0; +while(1) + d = (phi-1)*(xu - xl); + x1 = xl + d; + x2 = xu - d; + if f(x1,varargin{:}) < f(x2,varargin{:}) + xopt = x1; + xl = x2; + else + xopt = x2; + xu = x1; + end + iter=iter+1; + if xopt~=0, ea = (2 - phi) * abs((xu - xl) / xopt) * 100;end + if ea <= es | iter >= maxit,break,end +end +x=xopt;fx=f(xopt,varargin{:}); diff --git a/lecture_08/lecture_08.ipynb b/lecture_08/lecture_08.ipynb index 134ed8a..32b973c 100644 --- a/lecture_08/lecture_08.ipynb +++ b/lecture_08/lecture_08.ipynb @@ -2,41 +2,57 @@ "cells": [ { "cell_type": "code", - "execution_count": 2, + "execution_count": 6, "metadata": { - "collapsed": false + "collapsed": true }, "outputs": [], "source": [ - "%plot --format svg" + "function h = parabola(x)\n", + " y=1;\n", + " h=sum(x.^2+y.^2-1);\n", + "end\n", + " " ] }, { - "cell_type": "markdown", - "metadata": {}, + "cell_type": "code", + "execution_count": 9, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "ans = 0\r\n" + ] + } + ], "source": [ - "# Optimization\n", - "\n", - "Many problems involve finding a minimum or maximum based on given constraints. Engineers and scientists typically use energy balance equations to find the conditions of minimum energy, but this value may never go to 0 and the actual value of energy may not be of interest. " + "fzero(@(x) parabola(x),1)" ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 10, "metadata": { - "collapsed": true + "collapsed": false }, "outputs": [], - "source": [] + "source": [ + "%plot --format svg" + ] }, { - "cell_type": "code", - "execution_count": null, - "metadata": { - "collapsed": true - }, - "outputs": [], - "source": [] + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Optimization\n", + "\n", + "Many problems involve finding a minimum or maximum based on given constraints. Engineers and scientists typically use energy balance equations to find the conditions of minimum energy, but this value may never go to 0 and the actual value of energy may not be of interest. " + ] }, { "cell_type": "markdown", @@ -53,7 +69,7 @@ "Computational Tools to Study and Predict the Long-Term Stability of Nanowires.\n", "By Martin E. Zoloff Michoff, Patricio Vélez, Sergio A. Dassie and Ezequiel P. M. Leiva \n", "\n", - "![Model of Gold chain, from molecular dynamics simulation](Auchain_model.gif)\n", + "![Model of Gold chain, from molecular dynamics simulation](Auchain_model.png)\n", "\n", "[Single atom gold chain mechanics](http://www.uam.es/personal_pdi/ciencias/agrait/)\n", "\n", @@ -64,7 +80,7 @@ }, { "cell_type": "code", - "execution_count": 3, + "execution_count": 11, "metadata": { "collapsed": false }, @@ -243,7 +259,7 @@ }, { "cell_type": "code", - "execution_count": 4, + "execution_count": 12, "metadata": { "collapsed": false }, @@ -281,7 +297,7 @@ }, { "cell_type": "code", - "execution_count": 5, + "execution_count": 13, "metadata": { "collapsed": false }, @@ -515,7 +531,7 @@ }, { "cell_type": "code", - "execution_count": 6, + "execution_count": 14, "metadata": { "collapsed": false }, @@ -738,7 +754,7 @@ }, { "cell_type": "code", - "execution_count": 10, + "execution_count": 15, "metadata": { "collapsed": false }, @@ -748,7 +764,7 @@ "output_type": "stream", "text": [ "current_min = -0.038904\n", - "error_app = 0\n" + "error_app = 0.13359\n" ] }, { @@ -888,30 +904,30 @@ "\tgnuplot_plot_2a\n", "\n", "\t \n", - "\t\n", + "\t\n", "\n", "\t\n", "\tgnuplot_plot_3a\n", "\n", "\t \n", - "\t\n", + "\t\n", "\n", "\t\n", "\tgnuplot_plot_4a\n", "\n", "\t \n", - "\t\n", + "\t\n", "\n", "\t\n", "\tgnuplot_plot_5a\n", "\n", "\t \n", - "\t\n", + "\t\n", "\n", "\t\n", "\tgnuplot_plot_6a\n", "\n", - "\t\n", + "\t\n", "\t\n", "\n", "\n", @@ -961,7 +977,7 @@ }, { "cell_type": "code", - "execution_count": 13, + "execution_count": 16, "metadata": { "collapsed": false }, @@ -970,7 +986,7 @@ "name": "stdout", "output_type": "stream", "text": [ - "current_min = -0.039000\n", + "current_min = -0.038904\n", "error_app = 0\n" ] }, @@ -1111,30 +1127,30 @@ "\tgnuplot_plot_2a\n", "\n", "\t \n", - "\t\n", + "\t\n", "\n", "\t\n", "\tgnuplot_plot_3a\n", "\n", "\t \n", - "\t\n", + "\t\n", "\n", "\t\n", "\tgnuplot_plot_4a\n", "\n", "\t \n", - "\t\n", + "\t\n", "\n", "\t\n", "\tgnuplot_plot_5a\n", "\n", "\t \n", - "\t\n", + "\t\n", "\n", "\t\n", "\tgnuplot_plot_6a\n", "\n", - "\t\n", + "\t\n", "\t\n", "\n", "\n", @@ -1188,21 +1204,18 @@ "collapsed": true }, "source": [ - "## Parabolc Interpolation\n", + "## Parabolic Interpolation\n", "\n", "Near a minimum/maximum, the function resembles a parabola. With three data points, it is possible to fit a parabola to the function. So, given a lower and upper bound, we can choose the midpoint and fit a parabola to the 3 x,f(x) coordinates. \n", "\n", - "$ x_{4} =x_{2} - \\frac{1}{2} \\frac{ \\left(x_{2} -x_{1} \\right)^{2} \\left[f\n", - " \\left(x_{2} \\right)-f \\left(x_{3} \\right) \\right]- \\left(x_{2} -x_{3} \\right)^{2}\n", - " \\left[f \\left(x_{2} \\right)-f \\left(x_{1} \\right) \\right]}{ \\left(x_{2} -x_{1} \\right) \\left[f \\left(x_{2} \\right)-f \\left(x_{3} \\right)\n", - " \\right]- \\left(x_{2} -x_{3} \\right) \\left[f \\left(x_{2} \\right)-f \\left(x_{1} \\right) \\right]}$\n", + "$$ x_{4} =x_{2} - \\frac{1}{2}\\frac{\\left(x_{2} -x_{1}\\right)^{2} \\left(f\\left(x_{2} \\right)-f \\left(x_{3} \\right) \\right)- \\left(x_{2} -x_{3} \\right)^{2}\\left(f \\left(x_{2} \\right)-f \\left(x_{1} \\right) \\right)}{ \\left(x_{2} -x_{1} \\right) \\left(f \\left(x_{2} \\right)-f \\left(x_{3} \\right)\\right)- \\left(x_{2} -x_{3} \\right) \\left(f \\left(x_{2} \\right)-f \\left(x_{1} \\right) \\right)}$$\n", "\n", "Where $x_{4}$ location of the maximum or minimum of the parabola. " ] }, { "cell_type": "code", - "execution_count": 31, + "execution_count": 23, "metadata": { "collapsed": false }, @@ -1436,18 +1449,11 @@ }, { "cell_type": "code", - "execution_count": 32, + "execution_count": 24, "metadata": { "collapsed": false }, "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "f3 = -0.038878\r\n" - ] - }, { "data": { "image/svg+xml": [ @@ -1653,18 +1659,11 @@ }, { "cell_type": "code", - "execution_count": 33, + "execution_count": 25, "metadata": { "collapsed": false }, "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "f3 = -0.038878\r\n" - ] - }, { "data": { "image/svg+xml": [ @@ -1883,12 +1882,14 @@ "\n", "Now, we have two routines to find minimums of a univariate function (Golden Ratio and Parabolic). Let's use these to solve for the minimum energy given a range of applied forces to the single atom gold chain\n", "\n", - "$E_{total}(\\Delta x) = E_{LJ}(x_{min}+\\Delta x) - F \\cdot \\Delta x$" + "$E_{total}(\\Delta x) = E_{LJ}(x_{min}+\\Delta x) - F \\cdot \\Delta x$\n", + "\n", + "$1 aJ = 10^{-18} J = 1~nN* 1~nm = 10^{-9}N * 10^{-9} m$" ] }, { "cell_type": "code", - "execution_count": 37, + "execution_count": 26, "metadata": { "collapsed": false }, @@ -1897,14 +1898,834 @@ "name": "stdout", "output_type": "stream", "text": [ - "xmin = 3.2933\r\n" + "warning: Matlab-style short-circuit operation performed for operator |\n", + "warning: called from\n", + " goldmin at line 17 column 1\n", + "warning: Matlab-style short-circuit operation performed for operator |\n", + "warning: Matlab-style short-circuit operation performed for operator |\n", + "xmin = 0.32933\n", + "Emin = -2.7096e-04\n" ] + }, + { + "data": { + "image/svg+xml": [ + "\n", + "\n", + "Gnuplot\n", + "Produced by GNUPLOT 5.0 patchlevel 3 \n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t \n", + "\t \n", + "\t\n", + "\t\n", + "\t \n", + "\t \n", + "\t\n", + "\n", + "\n", + "\n", + "\n", + "\t\n", + "\t\t\n", + "\t\n", + "\n", + "\n", + "\n", + "\n", + "\t\t\n", + "\t\t-0.0004\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t-0.0002\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.0002\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.0004\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.0006\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.2\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.3\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.4\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.5\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.6\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.7\n", + "\t\n", + "\n", + "\n", + "\n", + "\n", + "\t\n", + "\n", + "\t\n", + "\t\tLennard Jones Potential (aJ/atom)\n", + "\t\n", + "\n", + "\n", + "\t\n", + "\t\tbond length (nm)\n", + "\t\n", + "\n", + "\n", + "\n", + "\tgnuplot_plot_1a\n", + "\n", + "\n", + "\n", + "\t\n", + "\t\n", + "\tgnuplot_plot_2a\n", + "\n", + "\t \n", + "\t\n", + "\n", + "\t\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "" + ], + "text/plain": [ + "" + ] + }, + "metadata": {}, + "output_type": "display_data" } ], "source": [ - "xmin = fminbnd(Au_x,2.8,6)\n", + "epsilon = 0.039; % kcal/mol\n", + "epsilon = epsilon*6.9477e-21; % J/atom\n", + "epsilon = epsilon*1e18; % aJ/J\n", + "% final units for epsilon are aJ\n", + "\n", + "sigma = 2.934; % Angstrom\n", + "sigma = sigma*0.10; % nm/Angstrom\n", + "x=linspace(2.8,6,200)*0.10; % bond length in um\n", "\n", - "Etotal = @(dx,F) Au_x(xmin+dx)-F.*dx;\n" + "Ex = lennard_jones(x,sigma,epsilon);\n", + "\n", + "%[Emin,imin]=min(Ex);\n", + "\n", + "[xmin,Emin] = goldmin(@(x) lennard_jones(x,sigma,epsilon),0.28,0.6)\n", + "\n", + "plot(x,Ex,xmin,Emin,'o')\n", + "ylabel('Lennard Jones Potential (aJ/atom)')\n", + "xlabel('bond length (nm)')\n", + "\n", + "Etotal = @(dx,F) lennard_jones(xmin+dx,sigma,epsilon)-F.*dx;\n" + ] + }, + { + "cell_type": "code", + "execution_count": 31, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "data": { + "image/svg+xml": [ + "\n", + "\n", + "Gnuplot\n", + "Produced by GNUPLOT 5.0 patchlevel 3 \n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t \n", + "\t \n", + "\t\n", + "\t\n", + "\t \n", + "\t \n", + "\t\n", + "\n", + "\n", + "\n", + "\n", + "\t\n", + "\t\t\n", + "\t\n", + "\n", + "\n", + "\n", + "\n", + "\t\t\n", + "\t\t0\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.0005\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.001\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.0015\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.002\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.0025\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.005\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.01\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.015\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.02\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.025\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.03\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.035\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.04\n", + "\t\n", + "\n", + "\n", + "\n", + "\n", + "\t\n", + "\n", + "\t\n", + "\t\tForce (nN)\n", + "\t\n", + "\n", + "\n", + "\t\n", + "\t\tdx (nm)\n", + "\t\n", + "\n", + "\n", + "\n", + "\tgnuplot_plot_1a\n", + "\n", + "\n", + "\n", + "\t\n", + "\t\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "" + ], + "text/plain": [ + "" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "N=50;\n", + "warning('off')\n", + "dx = zeros(1,N); % [in nm]\n", + "F_applied=linspace(0,0.0022,N); % [in nN]\n", + "for i=1:N\n", + " optmin=goldmin(@(dx) Etotal(dx,F_applied(i)),-0.001,0.035);\n", + " dx(i)=optmin;\n", + "end\n", + "\n", + "plot(dx,F_applied)\n", + "xlabel('dx (nm)')\n", + "ylabel('Force (nN)')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "For this function, it is possible to take a derivative and compare the analytical result:\n", + "\n", + "dE/dx = 0 = d(E_LJ)/dx - F\n", + "\n", + "F= d(E_LJ)/dx" + ] + }, + { + "cell_type": "code", + "execution_count": 32, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "F =\n", + "\n", + "@(dx) 4 * epsilon * 6 * (sigma ^ 6 ./ (xmin + dx) .^ 7 - 2 * sigma ^ 12 ./ (xmin + dx) .^ 13)\n", + "\n" + ] + }, + { + "data": { + "image/svg+xml": [ + "\n", + "\n", + "Gnuplot\n", + "Produced by GNUPLOT 5.0 patchlevel 3 \n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t \n", + "\t \n", + "\t\n", + "\t\n", + "\t \n", + "\t \n", + "\t\n", + "\n", + "\n", + "\n", + "\n", + "\t\n", + "\t\t\n", + "\t\n", + "\n", + "\n", + "\n", + "\n", + "\t\t\n", + "\t\t-0.0005\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.0005\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.001\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.0015\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.002\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.0025\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.01\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.02\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.03\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.04\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.05\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.06\n", + "\t\n", + "\n", + "\n", + "\n", + "\n", + "\t\n", + "\n", + "\n", + "\tgnuplot_plot_1a\n", + "\n", + "\n", + "\n", + "\t\n", + "\t\n", + "\tgnuplot_plot_2a\n", + "\n", + "\t\n", + "\t\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "" + ], + "text/plain": [ + "" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "dx_full=linspace(0,0.06,50);\n", + "F= @(dx) 4*epsilon*6*(sigma^6./(xmin+dx).^7-2*sigma^12./(xmin+dx).^13)\n", + "plot(dx_full,F(dx_full),dx,F_applied)" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "collapsed": true + }, + "source": [ + "## Curve-fitting\n", + "Another example is minimizing error in your approximation of a function. If you have data (now we have Force-displacement data) we can fit this to a function, such as:\n", + "\n", + "$F(x) = K_{1}\\Delta x + \\frac{1}{2} K_{2}(\\Delta x)^{2}$\n" + ] + }, + { + "cell_type": "code", + "execution_count": 33, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "function SSE = sse_of_parabola(K,xdata,ydata)\n", + " % calculate the sum of squares error for a parabola given a function, func, and xdata and ydata\n", + " % output is SSE=sum of squares error\n", + " K1=K(1);\n", + " K2=K(2);\n", + " y_function = K1*xdata+1/2*K2*xdata.^2;\n", + " SSE = sum((ydata-y_function).^2);\n", + "end\n", + " " + ] + }, + { + "cell_type": "code", + "execution_count": 34, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "\n", + "The nonlinear spring constants are K1=0.16 nN/nm and K2=-5.98 nN/nm^2\n", + "The mininum sum of squares error = 7.35e-08\n" + ] + } + ], + "source": [ + "[K,SSE_min]=fminsearch(@(K) sse_of_parabola(K,dx,F_applied),[1,1]);\n", + "fprintf('\\nThe nonlinear spring constants are K1=%1.2f nN/nm and K2=%1.2f nN/nm^2\\n',K)\n", + "fprintf('The mininum sum of squares error = %1.2e',SSE_min)" + ] + }, + { + "cell_type": "code", + "execution_count": 35, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "data": { + "image/svg+xml": [ + "\n", + "\n", + "Gnuplot\n", + "Produced by GNUPLOT 5.0 patchlevel 3 \n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t \n", + "\t \n", + "\t\n", + "\t\n", + "\t \n", + "\t \n", + "\t\n", + "\n", + "\n", + "\n", + "\n", + "\t\n", + "\t\t\n", + "\t\n", + "\n", + "\n", + "\n", + "\n", + "\t\t\n", + "\t\t0\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.0005\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.001\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.0015\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.002\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.0025\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.005\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.01\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.015\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.02\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.025\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.03\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.035\n", + "\t\n", + "\n", + "\n", + "\t\t\n", + "\t\t0.04\n", + "\t\n", + "\n", + "\n", + "\n", + "\n", + "\t\n", + "\n", + "\n", + "\tgnuplot_plot_1a\n", + "\n", + "\n", + "\n", + "\t \n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\t\n", + "\n", + "\t\n", + "\tgnuplot_plot_2a\n", + "\n", + "\t\n", + "\t\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "" + ], + "text/plain": [ + "" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "plot(dx,F_applied,'o',dx,K(1)*dx+1/2*K(2)*dx.^2)" ] }, { diff --git a/lecture_08/lecture_08.md b/lecture_08/lecture_08.md new file mode 100644 index 0000000..79edc0a --- /dev/null +++ b/lecture_08/lecture_08.md @@ -0,0 +1,626 @@ + + +```octave +%plot --format svg +``` + +# Optimization + +Many problems involve finding a minimum or maximum based on given constraints. Engineers and scientists typically use energy balance equations to find the conditions of minimum energy, but this value may never go to 0 and the actual value of energy may not be of interest. + +The Lennard-Jones potential is commonly used to model interatomic bonding. + +$E_{LJ}(x)=4\epsilon \left(\left(\frac{\sigma}{x}\right)^{12}-\left(\frac{\sigma}{x}\right)^{6}\right)$ + +Considering a 1-D gold chain, we can calculate the bond length, $x_{b}$, with no force applied to the chain and even for tension, F. This will allow us to calculate the nonlinear spring constant of a 1-D gold chain. + +![TEM image of Gold chain](au_chain.jpg) + +Computational Tools to Study and Predict the Long-Term Stability of Nanowires. +By Martin E. Zoloff Michoff, Patricio Vélez, Sergio A. Dassie and Ezequiel P. M. Leiva + +![Model of Gold chain, from molecular dynamics simulation](Auchain_model.png) + +[Single atom gold chain mechanics](http://www.uam.es/personal_pdi/ciencias/agrait/) + +### First, let's find the minimum energy $\min(E_{LJ}(x))$ + +## Brute force + + +```octave +setdefaults +epsilon = 0.039; % kcal/mol +sigma = 2.934; % Angstrom +x=linspace(2.8,6,200); % bond length in Angstrom + +Ex = lennard_jones(x,sigma,epsilon); + +[Emin,imin]=min(Ex); + +plot(x,Ex,x(imin),Emin,'o') +``` + + +![svg](lecture_08_files/lecture_08_3_0.svg) + + + +```octave +x(imin-1) +x(imin) +x(imin+1) +``` + + ans = 3.2824 + ans = 3.2985 + ans = 3.3146 + + +## Golden Search Algorithm + +We can't just look for a sign change for the problem (unless we can take a derivative) so we need a new approach to determine whether we have a maximum between the two bounds. + +Rather than using the midpoint of initial bounds, here the problem is more difficult. We need to compare the values of 4 function evaluations. The golden search uses the golden ratio to determine two interior points. + +![golden ratio](goldenratio.png) + +Start with bounds of 2.5 and 6 Angstrom. + + +```octave +% define Au atomic potential +epsilon = 0.039; % kcal/mol +sigma = 2.934; % Angstrom +Au_x= @(x) lennard_jones(x,sigma,epsilon); + +% calculate golden ratio +phi = 1/2+sqrt(5)/2; +% set initial limits +x_l=2.8; +x_u=6; + +% Iteration #1 +d=(phi-1)*(x_u-x_l); + +x1=x_l+d; % define point 1 +x2=x_u-d; % define point 2 + + +% evaluate Au_x(x1) and Au_x(x2) + +f1=Au_x(x1); +f2=Au_x(x2); +plot(x,Au_x(x),x_l,Au_x(x_l),'ro',x2,f2,'rs',x1,f1,'gs',x_u,Au_x(x_u),'go') +hold on; + +if f2x2 + plot(x4,f4,'*',[x1,x2],[f1,f2]) + x1=x2; + f1=f2; + else + plot(x4,f4,'*',[x3,x2],[f3,f2]) + x3=x2; + f3=f2; + end + x2=x4; f2=f4; +else + error('no minimum in bracket') +end +hold off +``` + + +![svg](lecture_08_files/lecture_08_11_0.svg) + + + +```octave +p = polyfit([x1,x2,x3],[f1,f2,f3],2); +x_fit = linspace(x1,x3,20); +y_fit = polyval(p,x_fit); + +plot(x,Au_x(x),x_fit,y_fit,[x1,x2,x3],[f1,f2,f3],'o') +hold on +if f2x2 + plot(x4,f4,'*',[x1,x2],[f1,f2]) + x1=x2; + f1=f2; + else + plot(x4,f4,'*',[x3,x2],[f3,f2]) + x3=x2; + f3=f2; + end + x2=x4; f2=f4; +else + error('no minimum in bracket') +end +hold off + +``` + + +![svg](lecture_08_files/lecture_08_12_0.svg) + + + +```octave +p = polyfit([x1,x2,x3],[f1,f2,f3],2); +x_fit = linspace(x1,x3,20); +y_fit = polyval(p,x_fit); + +plot(x,Au_x(x),x_fit,y_fit,[x1,x2,x3],[f1,f2,f3],'o') +hold on +if f2x2 + plot(x4,f4,'*',[x1,x2],[f1,f2]) + x1=x2; + f1=f2; + else + plot(x4,f4,'*',[x3,x2],[f3,f2]) + x3=x2; + f3=f2; + end + x2=x4; f2=f4; +else + error('no minimum in bracket') +end +hold off + +``` + + +![svg](lecture_08_files/lecture_08_13_0.svg) + + +Parabolic interpolation does not converge in many scenarios even though it it a bracketing method. Instead, functions like `fminbnd` in Matlab and Octave use a combination of the two (Golden Ratio and Parabolic) + +## Using the solutions to minimization for the nonlinear spring constant + +Now, we have two routines to find minimums of a univariate function (Golden Ratio and Parabolic). Let's use these to solve for the minimum energy given a range of applied forces to the single atom gold chain + +$E_{total}(\Delta x) = E_{LJ}(x_{min}+\Delta x) - F \cdot \Delta x$ + + +```octave +epsilon = 0.039; % kcal/mol +epsilon = epsilon*6.9477e-21; % J/atom +epsilon = epsilon*1e18; % aJ/J +% final units for epsilon are aJ + +sigma = 2.934; % Angstrom +sigma = sigma*0.10; % nm/Angstrom +x=linspace(2.8,6,200)*0.10; % bond length in um + +Ex = lennard_jones(x,sigma,epsilon); + +%[Emin,imin]=min(Ex); + +[xmin,Emin] = fminbnd(@(x) lennard_jones(x,sigma,epsilon),0.28,0.6) + +plot(x,Ex,xmin,Emin,'o') +ylabel('Lennard Jones Potential (aJ/atom)') +xlabel('bond length (nm)') + +Etotal = @(dx,F) lennard_jones(xmin+dx,sigma,epsilon)-F.*dx; + +``` + + xmin = 0.32933 + Emin = -2.7096e-04 + + + +![svg](lecture_08_files/lecture_08_16_1.svg) + + + +```octave +N=50; +dx = zeros(1,N); % [in nm] +F_applied=linspace(0,0.0022,N); % [in nN] +for i=1:N + optmin=goldmin(@(dx) Etotal(dx,F_applied(i)),-0.001,0.035); + dx(i)=optmin; +end + +plot(dx,F_applied) +xlabel('dx (nm)') +ylabel('Force (nN)') +``` + + warning: Matlab-style short-circuit operation performed for operator | + warning: called from + goldmin at line 17 column 1 + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + warning: Matlab-style short-circuit operation performed for operator | + + + +![svg](lecture_08_files/lecture_08_17_1.svg) + + +For this function, it is possible to take a derivative and compare the analytical result: + + +```octave +dx_full=linspace(0,0.06,50); +F= @(dx) 4*epsilon*6*(sigma^6./(xmin+dx).^7-2*sigma^12./(xmin+dx).^13) +plot(dx_full,F(dx_full),dx,F_applied) +``` + + F = + + @(dx) 4 * epsilon * 6 * (sigma ^ 6 ./ (xmin + dx) .^ 7 - 2 * sigma ^ 12 ./ (xmin + dx) .^ 13) + + + + +![svg](lecture_08_files/lecture_08_19_1.svg) + + +## Curve-fitting +Another example is minimizing error in your approximation of a function. If you have data (now we have Force-displacement data) we can fit this to a function, such as: + +$F(x) = K_{1}\Delta x + \frac{1}{2} K_{2}(\Delta x)^{2}$ + + + +```octave +function SSE = sse_of_parabola(K,xdata,ydata) + % calculate the sum of squares error for a parabola given a function, func, and xdata and ydata + % output is SSE=sum of squares error + K1=K(1); + K2=K(2); + y_function = K1*xdata+1/2*K2*xdata.^2; + SSE = sum((ydata-y_function).^2); +end + +``` + + +```octave +[K,SSE_min]=fminsearch(@(K) sse_of_parabola(K,dx,F_applied),[1,1]); +fprintf('\nThe nonlinear spring constants are K1=%1.2f nN/nm and K2=%1.2f nN/nm^2\n',K) +fprintf('The mininum sum of squares error = %1.2e',SSE_min) +``` + + + The nonlinear spring constants are K1=0.16 nN/nm and K2=-5.98 nN/nm^2 + The mininum sum of squares error = 7.35e-08 + + + +```octave +plot(dx,F_applied,'o',dx,K(1)*dx+1/2*K(2)*dx.^2) +``` + + +![svg](lecture_08_files/lecture_08_23_0.svg) + diff --git a/lecture_08/lecture_08.pdf b/lecture_08/lecture_08.pdf new file mode 100644 index 0000000..c251694 Binary files /dev/null and b/lecture_08/lecture_08.pdf differ diff --git a/lecture_08/lecture_08_files/lecture_08_11_0.svg b/lecture_08/lecture_08_files/lecture_08_11_0.svg new file mode 100644 index 0000000..c61f790 --- /dev/null +++ b/lecture_08/lecture_08_files/lecture_08_11_0.svg @@ -0,0 +1,168 @@ + + +Gnuplot +Produced by GNUPLOT 5.0 patchlevel 3 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + -0.06 + + + + + -0.04 + + + + + -0.02 + + + + + 0 + + + + + 0.02 + + + + + 0.04 + + + + + 0.06 + + + + + 0.08 + + + + + 2.5 + + + + + 3 + + + + + 3.5 + + + + + 4 + + + + + 4.5 + + + + + 5 + + + + + 5.5 + + + + + 6 + + + + + + + + + gnuplot_plot_1a + + + + + + gnuplot_plot_2a + + + + gnuplot_plot_3a + + + + + + + + gnuplot_plot_4a + + + + + + gnuplot_plot_5a + + + + + + + + + + + + + \ No newline at end of file diff --git a/lecture_08/lecture_08_files/lecture_08_12_0.svg b/lecture_08/lecture_08_files/lecture_08_12_0.svg new file mode 100644 index 0000000..b641032 --- /dev/null +++ b/lecture_08/lecture_08_files/lecture_08_12_0.svg @@ -0,0 +1,163 @@ + + +Gnuplot +Produced by GNUPLOT 5.0 patchlevel 3 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + -0.04 + + + + + -0.02 + + + + + 0 + + + + + 0.02 + + + + + 0.04 + + + + + 0.06 + + + + + 0.08 + + + + + 2.5 + + + + + 3 + + + + + 3.5 + + + + + 4 + + + + + 4.5 + + + + + 5 + + + + + 5.5 + + + + + 6 + + + + + + + + + gnuplot_plot_1a + + + + + + gnuplot_plot_2a + + + + gnuplot_plot_3a + + + + + + + + gnuplot_plot_4a + + + + + + gnuplot_plot_5a + + + + + + + + + + + + + \ No newline at end of file diff --git a/lecture_08/lecture_08_files/lecture_08_13_0.svg b/lecture_08/lecture_08_files/lecture_08_13_0.svg new file mode 100644 index 0000000..3e467e5 --- /dev/null +++ b/lecture_08/lecture_08_files/lecture_08_13_0.svg @@ -0,0 +1,163 @@ + + +Gnuplot +Produced by GNUPLOT 5.0 patchlevel 3 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + -0.04 + + + + + -0.02 + + + + + 0 + + + + + 0.02 + + + + + 0.04 + + + + + 0.06 + + + + + 0.08 + + + + + 2.5 + + + + + 3 + + + + + 3.5 + + + + + 4 + + + + + 4.5 + + + + + 5 + + + + + 5.5 + + + + + 6 + + + + + + + + + gnuplot_plot_1a + + + + + + gnuplot_plot_2a + + + + gnuplot_plot_3a + + + + + + + + gnuplot_plot_4a + + + + + + gnuplot_plot_5a + + + + + + + + + + + + + \ No newline at end of file diff --git a/lecture_08/lecture_08_files/lecture_08_16_1.svg b/lecture_08/lecture_08_files/lecture_08_16_1.svg new file mode 100644 index 0000000..ff27c29 --- /dev/null +++ b/lecture_08/lecture_08_files/lecture_08_16_1.svg @@ -0,0 +1,142 @@ + + +Gnuplot +Produced by GNUPLOT 5.0 patchlevel 3 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + -0.0004 + + + + + -0.0002 + + + + + 0 + + + + + 0.0002 + + + + + 0.0004 + + + + + 0.0006 + + + + + 0.2 + + + + + 0.3 + + + + + 0.4 + + + + + 0.5 + + + + + 0.6 + + + + + 0.7 + + + + + + + + + Lennard Jones Potential (aJ/atom) + + + + + bond length (nm) + + + + + gnuplot_plot_1a + + + + + + gnuplot_plot_2a + + + + + + + + + + + + + + + \ No newline at end of file diff --git a/lecture_08/lecture_08_files/lecture_08_17_1.svg b/lecture_08/lecture_08_files/lecture_08_17_1.svg new file mode 100644 index 0000000..1f3103b --- /dev/null +++ b/lecture_08/lecture_08_files/lecture_08_17_1.svg @@ -0,0 +1,136 @@ + + +Gnuplot +Produced by GNUPLOT 5.0 patchlevel 3 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + 0 + + + + + 0.0005 + + + + + 0.001 + + + + + 0.0015 + + + + + 0.002 + + + + + 0.0025 + + + + + -0.01 + + + + + 0 + + + + + 0.01 + + + + + 0.02 + + + + + 0.03 + + + + + 0.04 + + + + + + + + + Force (nN) + + + + + dx (nm) + + + + + gnuplot_plot_1a + + + + + + + + + + + + + + + \ No newline at end of file diff --git a/lecture_08/lecture_08_files/lecture_08_19_1.svg b/lecture_08/lecture_08_files/lecture_08_19_1.svg new file mode 100644 index 0000000..2a891df --- /dev/null +++ b/lecture_08/lecture_08_files/lecture_08_19_1.svg @@ -0,0 +1,135 @@ + + +Gnuplot +Produced by GNUPLOT 5.0 patchlevel 3 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + 0 + + + + + 0.0005 + + + + + 0.001 + + + + + 0.0015 + + + + + 0.002 + + + + + 0.0025 + + + + + 0 + + + + + 0.01 + + + + + 0.02 + + + + + 0.03 + + + + + 0.04 + + + + + 0.05 + + + + + 0.06 + + + + + + + + + gnuplot_plot_1a + + + + + + gnuplot_plot_2a + + + + + + + + + + + + + \ No newline at end of file diff --git a/lecture_08/lecture_08_files/lecture_08_23_0.svg b/lecture_08/lecture_08_files/lecture_08_23_0.svg new file mode 100644 index 0000000..8fa3bb1 --- /dev/null +++ b/lecture_08/lecture_08_files/lecture_08_23_0.svg @@ -0,0 +1,196 @@ + + +Gnuplot +Produced by GNUPLOT 5.0 patchlevel 3 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + 0 + + + + + 0.0005 + + + + + 0.001 + + + + + 0.0015 + + + + + 0.002 + + + + + 0.0025 + + + + + 0 + + + + + 0.005 + + + + + 0.01 + + + + + 0.015 + + + + + 0.02 + + + + + 0.025 + + + + + 0.03 + + + + + 0.035 + + + + + 0.04 + + + + + + + + + gnuplot_plot_1a + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + gnuplot_plot_2a + + + + + + + + + + + + + \ No newline at end of file diff --git a/lecture_08/lecture_08_files/lecture_08_3_0.svg b/lecture_08/lecture_08_files/lecture_08_3_0.svg new file mode 100644 index 0000000..2c38150 --- /dev/null +++ b/lecture_08/lecture_08_files/lecture_08_3_0.svg @@ -0,0 +1,147 @@ + + +Gnuplot +Produced by GNUPLOT 5.0 patchlevel 3 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + -0.04 + + + + + -0.02 + + + + + 0 + + + + + 0.02 + + + + + 0.04 + + + + + 0.06 + + + + + 0.08 + + + + + 2.5 + + + + + 3 + + + + + 3.5 + + + + + 4 + + + + + 4.5 + + + + + 5 + + + + + 5.5 + + + + + 6 + + + + + + + + + gnuplot_plot_1a + + + + + + gnuplot_plot_2a + + + + + + + + + + + + + + + \ No newline at end of file diff --git a/lecture_08/lecture_08_files/lecture_08_6_1.svg b/lecture_08/lecture_08_files/lecture_08_6_1.svg new file mode 100644 index 0000000..ecf7a18 --- /dev/null +++ b/lecture_08/lecture_08_files/lecture_08_6_1.svg @@ -0,0 +1,169 @@ + + +Gnuplot +Produced by GNUPLOT 5.0 patchlevel 3 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + -0.04 + + + + + -0.02 + + + + + 0 + + + + + 0.02 + + + + + 0.04 + + + + + 0.06 + + + + + 0.08 + + + + + 2.5 + + + + + 3 + + + + + 3.5 + + + + + 4 + + + + + 4.5 + + + + + 5 + + + + + 5.5 + + + + + 6 + + + + + + + + + gnuplot_plot_1a + + + + + + gnuplot_plot_2a + + + + + + gnuplot_plot_3a + + + + + + gnuplot_plot_4a + + + + + + gnuplot_plot_5a + + + + + + gnuplot_plot_6a + + + + + + + + + + + + + \ No newline at end of file diff --git a/lecture_08/lecture_08_files/lecture_08_7_1.svg b/lecture_08/lecture_08_files/lecture_08_7_1.svg new file mode 100644 index 0000000..e63011c --- /dev/null +++ b/lecture_08/lecture_08_files/lecture_08_7_1.svg @@ -0,0 +1,169 @@ + + +Gnuplot +Produced by GNUPLOT 5.0 patchlevel 3 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + -0.04 + + + + + -0.02 + + + + + 0 + + + + + 0.02 + + + + + 0.04 + + + + + 0.06 + + + + + 0.08 + + + + + 2.5 + + + + + 3 + + + + + 3.5 + + + + + 4 + + + + + 4.5 + + + + + 5 + + + + + 5.5 + + + + + 6 + + + + + + + + + gnuplot_plot_1a + + + + + + gnuplot_plot_2a + + + + + + gnuplot_plot_3a + + + + + + gnuplot_plot_4a + + + + + + gnuplot_plot_5a + + + + + + gnuplot_plot_6a + + + + + + + + + + + + + \ No newline at end of file diff --git a/lecture_08/lecture_08_files/lecture_08_8_1.svg b/lecture_08/lecture_08_files/lecture_08_8_1.svg new file mode 100644 index 0000000..873b04d --- /dev/null +++ b/lecture_08/lecture_08_files/lecture_08_8_1.svg @@ -0,0 +1,169 @@ + + +Gnuplot +Produced by GNUPLOT 5.0 patchlevel 3 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + -0.04 + + + + + -0.02 + + + + + 0 + + + + + 0.02 + + + + + 0.04 + + + + + 0.06 + + + + + 0.08 + + + + + 2.5 + + + + + 3 + + + + + 3.5 + + + + + 4 + + + + + 4.5 + + + + + 5 + + + + + 5.5 + + + + + 6 + + + + + + + + + gnuplot_plot_1a + + + + + + gnuplot_plot_2a + + + + + + gnuplot_plot_3a + + + + + + gnuplot_plot_4a + + + + + + gnuplot_plot_5a + + + + + + gnuplot_plot_6a + + + + + + + + + + + + + \ No newline at end of file diff --git a/lecture_08/lecture_08_files/lecture_08_9_1.svg b/lecture_08/lecture_08_files/lecture_08_9_1.svg new file mode 100644 index 0000000..28fd681 --- /dev/null +++ b/lecture_08/lecture_08_files/lecture_08_9_1.svg @@ -0,0 +1,169 @@ + + +Gnuplot +Produced by GNUPLOT 5.0 patchlevel 3 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + -0.04 + + + + + -0.02 + + + + + 0 + + + + + 0.02 + + + + + 0.04 + + + + + 0.06 + + + + + 0.08 + + + + + 2.5 + + + + + 3 + + + + + 3.5 + + + + + 4 + + + + + 4.5 + + + + + 5 + + + + + 5.5 + + + + + 6 + + + + + + + + + gnuplot_plot_1a + + + + + + gnuplot_plot_2a + + + + + + gnuplot_plot_3a + + + + + + gnuplot_plot_4a + + + + + + gnuplot_plot_5a + + + + + + gnuplot_plot_6a + + + + + + + + + + + + + \ No newline at end of file diff --git a/lecture_08/octave-workspace b/lecture_08/octave-workspace index f505b0a..8c437bb 100644 Binary files a/lecture_08/octave-workspace and b/lecture_08/octave-workspace differ diff --git a/lecture_09/.ipynb_checkpoints/lecture_09-checkpoint.ipynb b/lecture_09/.ipynb_checkpoints/lecture_09-checkpoint.ipynb new file mode 100644 index 0000000..2fd6442 --- /dev/null +++ b/lecture_09/.ipynb_checkpoints/lecture_09-checkpoint.ipynb @@ -0,0 +1,6 @@ +{ + "cells": [], + "metadata": {}, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/lecture_09/.lecture_09.md.swp b/lecture_09/.lecture_09.md.swp new file mode 100644 index 0000000..4000b4f Binary files /dev/null and b/lecture_09/.lecture_09.md.swp differ diff --git a/lecture_09/even_odd_perm.svg b/lecture_09/even_odd_perm.svg new file mode 100644 index 0000000..eeaca19 --- /dev/null +++ b/lecture_09/even_odd_perm.svg @@ -0,0 +1,331 @@ + + + + + + + + + + + + + + + + + + + + + + + image/svg+xml + + + + + + + 1 + 2 + 3 + 4 + 5 + 6 + + + + 1 + 2 + 3 + 4 + 5 + 6 + + + + + even+ + odd- + + + diff --git a/lecture_09/lecture_09.ipynb b/lecture_09/lecture_09.ipynb new file mode 100644 index 0000000..4889d52 --- /dev/null +++ b/lecture_09/lecture_09.ipynb @@ -0,0 +1,1151 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Linear Algebra (Review/Introduction)\n", + "\n", + "Representation of linear equations:\n", + "\n", + "1. $5x_{1}+3x_{2} =1$\n", + "\n", + "2. $x_{1}+2x_{2}+3x_{3} =2$\n", + "\n", + "3. $x_{1}+x_{2}+x_{3} =3$\n", + "\n", + "in matrix form:\n", + "\n", + "$\\left[ \\begin{array}{ccc}\n", + "5 & 3 & 0 \\\\\n", + "1 & 2 & 3 \\\\\n", + "1 & 1 & 1 \\end{array} \\right]\n", + "\\left[\\begin{array}{c} \n", + "x_{1} \\\\ \n", + "x_{2} \\\\\n", + "x_{3}\\end{array}\\right]=\\left[\\begin{array}{c} \n", + "1 \\\\\n", + "2 \\\\\n", + "3\\end{array}\\right]$\n", + "\n", + "$Ax=b$" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Vectors \n", + "\n", + "column vector x (length of 3):\n", + "\n", + "$\\left[\\begin{array}{c} \n", + "x_{1} \\\\ \n", + "x_{2} \\\\\n", + "x_{3}\\end{array}\\right]$\n", + "\n", + "row vector y (length of 3):\n", + "\n", + "$\\left[\\begin{array}{ccc} \n", + "y_{1}~ \n", + "y_{2}~ \n", + "y_{3}\\end{array}\\right]$\n", + "\n", + "vector of length N:\n", + "\n", + "$\\left[\\begin{array}{c} \n", + "x_{1} \\\\ \n", + "x_{2} \\\\\n", + "\\vdots \\\\\n", + "x_{N}\\end{array}\\right]$\n", + "\n", + "The $i^{th}$ element of x is $x_{i}$" + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "x =\n", + "\n", + " 1 2 3 4 5 6 7 8 9 10\n", + "\n" + ] + } + ], + "source": [ + "x=[1:10]" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "ans =\n", + "\n", + " 1\n", + " 2\n", + " 3\n", + " 4\n", + " 5\n", + " 6\n", + " 7\n", + " 8\n", + " 9\n", + " 10\n", + "\n" + ] + } + ], + "source": [ + "x'" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Matrices\n", + "\n", + "Matrix A is 3x3:\n", + "\n", + "$A=\\left[ \\begin{array}{ccc}\n", + "5 & 3 & 0 \\\\\n", + "1 & 2 & 3 \\\\\n", + "1 & 1 & 1 \\end{array} \\right]$\n", + "\n", + "elements in the matrix are denoted $A_{row~column}$, $A_{23}=3$\n", + "\n", + "In general, matrix, B, can be any size, $M \\times N$ (M-rows and N-columns):\n", + "\n", + "$B=\\left[ \\begin{array}{cccc}\n", + "B_{11} & B_{12} &...& B_{1N} \\\\\n", + "B_{21} & B_{22} &...& B_{2N} \\\\\n", + "\\vdots & \\vdots &\\ddots& \\vdots \\\\\n", + "B_{M1} & B_{M2} &...& B_{MN}\\end{array} \\right]$" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "A =\n", + "\n", + " 5 3 0\n", + " 1 2 3\n", + " 1 1 1\n", + "\n" + ] + } + ], + "source": [ + "A=[5,3,0;1,2,3;1,1,1]" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Multiplication\n", + "\n", + "A column vector is a $1\\times N$ matrix and a row vector is a $M\\times 1$ matrix\n", + "\n", + "**Multiplying matrices is not commutative**\n", + "\n", + "$A B \\neq B A$\n", + "\n", + "Inner dimensions must agree, output is outer dimensions. \n", + "\n", + "A is $M1 \\times N1$ and B is $M2 \\times N2$\n", + "\n", + "C=AB\n", + "\n", + "Therefore N1=M2 and C is $M1 \\times N2$\n", + "\n", + "If $C'=BA$, then N2=M1 and C is $M2 \\times N1$\n", + "\n", + "\n", + "e.g. \n", + "$A=\\left[ \\begin{array}{cc}\n", + "5 & 3 & 0 \\\\\n", + "1 & 2 & 3 \\end{array} \\right]$\n", + "\n", + "$B=\\left[ \\begin{array}{cccc}\n", + "1 & 2 & 3 & 4 \\\\\n", + "5 & 6 & 7 & 8 \\\\\n", + "9 & 10 & 11 & 12 \\end{array} \\right]$\n", + "\n", + "C=AB\n", + "\n", + "$[2\\times 4] = [2 \\times 3][3 \\times 4]$\n", + "\n", + "The rule for multiplying matrices, A and B, is\n", + "\n", + "$C_{ij} = \\sum_{k=1}^{n}A_{ik}B_{kj}$\n", + "\n", + "In the previous example, \n", + "\n", + "$C_{11} = A_{11}B_{11}+A_{12}B_{21}+A_{13}B_{31} = 5*1+3*5+0*9=20$\n", + "\n", + "\n", + "Multiplication is associative:\n", + "\n", + "$(AB)C = A(BC)$\n", + "\n", + "and distributive:\n", + "\n", + "$A(B+C)=AB+AC$" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "A =\n", + "\n", + " 5 3 0\n", + " 1 2 3\n", + "\n", + "B =\n", + "\n", + " 1 2 3 4\n", + " 5 6 7 8\n", + " 9 10 11 12\n", + "\n" + ] + } + ], + "source": [ + "A=[5,3,0;1,2,3] \n", + "B=[1,2,3,4;5,6,7,8;9,10,11,12]" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "C =\n", + "\n", + " 0 0 0 0\n", + " 0 0 0 0\n", + "\n", + "C =\n", + "\n", + " 20 28 36 44\n", + " 38 44 50 56\n", + "\n" + ] + } + ], + "source": [ + "C=zeros(2,4)\n", + "C(1,1)=A(1,1)*B(1,1)+A(1,2)*B(2,1)+A(1,3)*B(3,1);\n", + "C(1,2)=A(1,1)*B(1,2)+A(1,2)*B(2,2)+A(1,3)*B(3,2);\n", + "\n", + "C=A*B" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "error: operator *: nonconformant arguments (op1 is 3x4, op2 is 2x3)\r\n" + ] + } + ], + "source": [ + "Cp=B*A" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Representation of a problem with Matrices and Vectors\n", + "\n", + "If you have a set of known output, $y_{1},~y_{2},~...y_{N}$ and a set of equations that\n", + "relate unknown inputs, $x_{1},~x_{2},~...x_{N}$, then these can be written in a vector\n", + "matrix format as:\n", + "\n", + "$y=Ax=\\left[\\begin{array}{c} \n", + "y_{1} \\\\ \n", + "y_{2} \\\\\n", + "\\vdots \\\\\n", + "y_{N}\\end{array}\\right]\n", + "=\n", + "A\\left[\\begin{array}{c} \n", + "x_{1} \\\\ \n", + "x_{2} \\\\\n", + "\\vdots \\\\\n", + "x_{N}\\end{array}\\right]$\n", + "\n", + "$A=\n", + "\\left[\\begin{array}{cccc} \n", + "| & | & & | \\\\ \n", + "a_{1} & a_{2} & ... & a_{N} \\\\ \n", + "| & | & & | \\end{array}\\right]$\n", + "\n", + "or \n", + "\n", + "$y = a_{1}x_{1} + a_{2}x_{2} +...+a_{N}x_{N}$\n", + "\n", + "where each $a_{i}$ is a column vector and $x_{i}$ is a scalar taken from the $i^{th}$\n", + "element of x.\n", + "\n", + "Consider the following problem, 4 masses are connected in series to 4 springs with K=10 N/m. What are the final positions of the masses? \n", + "\n", + "![Springs-masses](mass_springs.svg)\n", + "\n", + "The masses haves the following amounts, 1, 2, 3, and 4 kg for masses 1-4. Using a FBD for each mass:\n", + "\n", + "$m_{1}g+k(x_{2}-x_{1})-kx_{1}=0$\n", + "\n", + "$m_{2}g+k(x_{3}-x_{2})-k(x_{2}-x_{1})=0$\n", + "\n", + "$m_{3}g+k(x_{4}-x_{3})-k(x_{3}-x_{2})=0$\n", + "\n", + "$m_{4}g-k(x_{4}-x_{3})=0$\n", + "\n", + "in matrix form:\n", + "\n", + "$\\left[ \\begin{array}{cccc}\n", + "2k & -k & 0 & 0 \\\\\n", + "-k & 2k & -k & 0 \\\\\n", + "0 & -k & 2k & -k \\\\\n", + "0 & 0 & -k & k \\end{array} \\right]\n", + "\\left[ \\begin{array}{c}\n", + "x_{1} \\\\\n", + "x_{2} \\\\\n", + "x_{3} \\\\\n", + "x_{4} \\end{array} \\right]=\n", + "\\left[ \\begin{array}{c}\n", + "m_{1}g \\\\\n", + "m_{2}g \\\\\n", + "m_{3}g \\\\\n", + "m_{4}g \\end{array} \\right]$" + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "K =\n", + "\n", + " 20 -10 0 0\n", + " -10 20 -10 0\n", + " 0 -10 20 -10\n", + " 0 0 -10 10\n", + "\n", + "y =\n", + "\n", + " 9.8100\n", + " 19.6200\n", + " 29.4300\n", + " 39.2400\n", + "\n", + "x =\n", + "\n", + " 9.8100\n", + " 18.6390\n", + " 25.5060\n", + " 29.4300\n", + "\n" + ] + } + ], + "source": [ + "k=10; % N/m\n", + "m1=1; % kg\n", + "m2=2;\n", + "m3=3;\n", + "m4=4;\n", + "g=9.81; % m/s^2\n", + "K=[2*k -k 0 0; -k 2*k -k 0; 0 -k 2*k -k; 0 0 -k k]\n", + "y=[m1*g;m2*g;m3*g;m4*g]\n", + "\n", + "x=K\\y" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Matrix Operations\n", + "\n", + "Identity matrix `eye(M,N)` **Assume M=N unless specfied**\n", + "\n", + "$I_{ij} =1$ if $i=j$ \n", + "\n", + "$I_{ij} =0$ if $i\\neq j$ \n", + "\n", + "AI=A=IA\n", + "\n", + "Diagonal matrix, D, is similar to the identity matrix, but each diagonal element can be any\n", + "number\n", + "\n", + "$D_{ij} =d_{i}$ if $i=j$ \n", + "\n", + "$D_{ij} =0$ if $i\\neq j$ " + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Transpose\n", + "\n", + "The transpose of a matrix changes the rows -> columns and columns-> rows\n", + "\n", + "$(A^{T})_{ij} = A_{ji}$\n", + "\n", + "for diagonal matrices, $D^{T}=D$\n", + "\n", + "in general, the transpose has the following qualities:\n", + " \n", + "1. $(A^{T})^{T}=A$\n", + "\n", + "2. $(AB)^{T} = B^{T}A^{T}$\n", + "\n", + "3. $(A+B)^{T} = A^{T} +B^{T}$ \n", + "\n", + "All matrices have a symmetric and antisymmetric part:\n", + "\n", + "$A= 1/2(A+A^{T}) +1/2(A-A^{T})$\n", + "\n", + "If $A=A^{T}$, then A is symmetric, if $A=-A^{T}$, then A is antisymmetric." + ] + }, + { + "cell_type": "code", + "execution_count": 18, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "ans =\n", + "\n", + " 20 38\n", + " 28 44\n", + " 36 50\n", + " 44 56\n", + "\n", + "ans =\n", + "\n", + " 20 38\n", + " 28 44\n", + " 36 50\n", + " 44 56\n", + "\n", + "error: operator *: nonconformant arguments (op1 is 3x2, op2 is 4x3)\n" + ] + } + ], + "source": [ + "(A*B)'\n", + "\n", + "B'*A' % if this wasnt true, then inner dimensions wouldn;t match\n", + "\n", + "A'*B' == (A*B)'" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Square matrix operations\n", + "\n", + "### Trace\n", + "\n", + "The trace of a square matrix is the sum of the diagonal elements\n", + "\n", + "$tr~A=\\sum_{i=1}^{N}A_{ii}$\n", + "\n", + "The trace is invariant, meaning that if you change the basis through a rotation, then the\n", + "trace remains the same. \n", + "\n", + "It also has the following properties:\n", + "\n", + "1. $tr~A=tr~A^{T}$\n", + "\n", + "2. $tr~A+tr~B=tr(A+B)$\n", + "\n", + "3. $tr(kA)=k tr~A$\n", + "\n", + "4. $tr(AB)=tr(BA)$" + ] + }, + { + "cell_type": "code", + "execution_count": 19, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "id_m =\n", + "\n", + "Diagonal Matrix\n", + "\n", + " 1 0 0\n", + " 0 1 0\n", + " 0 0 1\n", + "\n", + "ans = 3\n" + ] + } + ], + "source": [ + "id_m=eye(3)\n", + "trace(id_m)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Inverse\n", + "\n", + "The inverse of a square matrix, $A^{-1}$ is defined such that\n", + "\n", + "$A^{-1}A=I=AA^{-1}$\n", + "\n", + "Not all square matrices have an inverse, they can be *singular* or *non-invertible*\n", + "\n", + "The inverse has the following properties:\n", + "\n", + "1. $(A^{-1})^{-1}=A$\n", + "\n", + "2. $(AB)^{-1}=B^{-1}A^{-1}$\n", + "\n", + "3. $(A^{-1})^{T}=(A^{T})^{-1}$" + ] + }, + { + "cell_type": "code", + "execution_count": 20, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "A =\n", + "\n", + " 0.5762106 0.3533174 0.7172134\n", + " 0.7061664 0.4863733 0.9423436\n", + " 0.4255961 0.0016613 0.3561407\n", + "\n" + ] + } + ], + "source": [ + "A=rand(3,3)" + ] + }, + { + "cell_type": "code", + "execution_count": 22, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Ainv =\n", + "\n", + " 41.5613 -30.1783 -3.8467\n", + " 36.2130 -24.2201 -8.8415\n", + " -49.8356 36.1767 7.4460\n", + "\n" + ] + } + ], + "source": [ + "Ainv=inv(A)" + ] + }, + { + "cell_type": "code", + "execution_count": 26, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "B =\n", + "\n", + " 0.524529 0.470856 0.708116\n", + " 0.084491 0.730986 0.528292\n", + " 0.303545 0.782522 0.389282\n", + "\n" + ] + } + ], + "source": [ + "B=rand(3,3)" + ] + }, + { + "cell_type": "code", + "execution_count": 28, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "ans =\n", + "\n", + " -182.185 125.738 40.598\n", + " -133.512 97.116 17.079\n", + " 282.422 -200.333 -46.861\n", + "\n", + "ans =\n", + "\n", + " -182.185 125.738 40.598\n", + " -133.512 97.116 17.079\n", + " 282.422 -200.333 -46.861\n", + "\n", + "ans =\n", + "\n", + " 41.5613 36.2130 -49.8356\n", + " -30.1783 -24.2201 36.1767\n", + " -3.8467 -8.8415 7.4460\n", + "\n", + "ans =\n", + "\n", + " 41.5613 36.2130 -49.8356\n", + " -30.1783 -24.2201 36.1767\n", + " -3.8467 -8.8415 7.4460\n", + "\n" + ] + } + ], + "source": [ + "inv(A*B)\n", + "inv(B)*inv(A)\n", + "\n", + "inv(A')\n", + "\n", + "inv(A)'" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Orthogonal Matrices\n", + "\n", + "Vectors are *orthogonal* if $x^{T}$ y=0, and a vector is *normalized* if $||x||_{2}=1$. A\n", + "square matrix is *orthogonal* if all its column vectors are orthogonal to each other and\n", + "normalized. The column vectors are then called *orthonormal* and the following results\n", + "\n", + "$U^{T}U=I=UU^{T}$\n", + "\n", + "and \n", + "\n", + "$||Ux||_{2}=||x||_{2}$\n", + "\n", + "### Determinant\n", + "\n", + "The **determinant** of a matrix has 3 properties\n", + "\n", + "1. The determinant of the identity matrix is one, $|I|=1$\n", + "\n", + "2. If you multiply a single row by scalar $t$ then the determinant is $t|A|$:\n", + "\n", + "$t|A|=\\left[ \\begin{array}{cccc}\n", + "tA_{11} & tA_{12} &...& tA_{1N} \\\\\n", + "A_{21} & A_{22} &...& A_{2N} \\\\\n", + "\\vdots & \\vdots &\\ddots& \\vdots \\\\\n", + "A_{M1} & A_{M2} &...& A_{MN}\\end{array} \\right]$\n", + "\n", + "3. If you switch 2 rows, the determinant changes sign:\n", + "\n", + "\n", + "$-|A|=\\left[ \\begin{array}{cccc}\n", + "A_{21} & A_{22} &...& A_{2N} \\\\\n", + "A_{11} & A_{12} &...& A_{1N} \\\\\n", + "\\vdots & \\vdots &\\ddots& \\vdots \\\\\n", + "A_{M1} & A_{M2} &...& A_{MN}\\end{array} \\right]$\n", + "\n", + "4. inverse of the determinant is the determinant of the inverse:\n", + "\n", + "$|A^{-1}|=\\frac{1}{|A|}=|A|^{-1}$\n", + "\n", + "For a $2\\times2$ matrix, \n", + "\n", + "$|A|=\\left|\\left[ \\begin{array}{cc}\n", + "A_{11} & A_{12} \\\\\n", + "A_{21} & A_{22} \\\\\n", + "\\end{array} \\right]\\right| = A_{11}A_{22}-A_{21}A_{12}$\n", + "\n", + "For a $3\\times3$ matrix,\n", + "\n", + "$|A|=\\left|\\left[ \\begin{array}{ccc}\n", + "A_{11} & A_{12} & A_{13} \\\\\n", + "A_{21} & A_{22} & A_{23} \\\\\n", + "A_{31} & A_{32} & A_{33}\\end{array} \\right]\\right|=$\n", + "\n", + "$A_{11}A_{22}A_{33}+A_{12}A_{23}A_{31}+A_{13}A_{21}A_{32}\n", + "-A_{31}A_{22}A_{13}-A_{32}A_{23}A_{11}-A_{33}A_{21}A_{12}$\n", + "\n", + "For larger matrices, the determinant is \n", + "\n", + "$|A|=\n", + "\\frac{1}{n!}\n", + "\\sum_{i_{r}}^{N}\n", + "\\sum_{j_{r}}^{N}\n", + "\\epsilon_{i_{1}...i_{N}}\n", + "\\epsilon_{j_{1}...j_{N}}\n", + "A_{i_{1}j_{1}}\n", + "A_{i_{2}j_{2}}\n", + "...\n", + "A_{i_{N}j_{N}}$\n", + "\n", + "Where the Levi-Cevita symbol is $\\epsilon_{i_{1}i_{2}...i_{N}} = 1$ if it is an even\n", + "permutation and $\\epsilon_{i_{1}i_{2}...i_{N}} = -1$ if it is an odd permutation and\n", + "$\\epsilon_{i_{1}i_{2}...i_{N}} = 0$ if $i_{n}=i_{m}$ e.g. $i_{1}=i_{7}$\n", + "\n", + "![Levi-Cevita rule for even and odd permutations for an $N\\times N$\n", + "determinant](even_odd_perm.svg)\n", + "\n", + "So a $4\\times4$ matrix determinant,\n", + "\n", + "$|A|=\\left|\\left[ \\begin{array}{cccc}\n", + "A_{11} & A_{12} & A_{13} & A_{14} \\\\\n", + "A_{21} & A_{22} & A_{23} & A_{24} \\\\\n", + "A_{31} & A_{32} & A_{33} & A_{34} \\\\\n", + "A_{41} & A_{42} & A_{43} & A_{44} \\end{array} \\right]\\right|$\n", + "\n", + "$=\\epsilon_{1234}\\epsilon_{1234}A_{11}A_{22}A_{33}A_{44}+\n", + "\\epsilon_{2341}\\epsilon_{1234}A_{21}A_{32}A_{43}A_{14}+\n", + "\\epsilon_{3412}\\epsilon_{1234}A_{31}A_{42}A_{13}A_{24}+\n", + "\\epsilon_{4123}\\epsilon_{1234}A_{41}A_{12}A_{23}A_{34}+...\n", + "\\epsilon_{4321}\\epsilon_{4321}A_{44}A_{33}A_{22}A_{11}$\n", + "\n", + "### Special Case for determinants\n", + "\n", + "The determinant of a diagonal matrix $|D|=D_{11}D_{22}D_{33}...D_{NN}$. \n", + "\n", + "Similarly, if a matrix is upper triangular (so all values of $A_{ij}=0$ when $j +Babel <3.9q> and hyphenation patterns for 81 language(s) loaded. +! You can't use `macro parameter character #' in vertical mode. +l.2 # + Linear Algebra (Review/Introduction) +? +! Interruption. + + +l.2 # + Linear Algebra (Review/Introduction) +? xx + +Here is how much of TeX's memory you used: + 6 strings out of 493030 + 147 string characters out of 6136260 + 53067 words of memory out of 5000000 + 3641 multiletter control sequences out of 15000+600000 + 3640 words of font info for 14 fonts, out of 8000000 for 9000 + 1141 hyphenation exceptions out of 8191 + 5i,0n,1p,57b,8s stack positions out of 5000i,500n,10000p,200000b,80000s +No pages of output. diff --git a/lecture_09/lecture_09.md b/lecture_09/lecture_09.md new file mode 100644 index 0000000..a1d53d1 --- /dev/null +++ b/lecture_09/lecture_09.md @@ -0,0 +1,802 @@ + +# Linear Algebra (Review/Introduction) + +Representation of linear equations: + +1. $5x_{1}+3x_{2} =1$ + +2. $x_{1}+2x_{2}+3x_{3} =2$ + +3. $x_{1}+x_{2}+x_{3} =3$ + +in matrix form: + +$\left[ \begin{array}{ccc} +5 & 3 & 0 \\ +1 & 2 & 3 \\ +1 & 1 & 1 \end{array} \right] +\left[\begin{array}{c} +x_{1} \\ +x_{2} \\ +x_{3}\end{array}\right]=\left[\begin{array}{c} +1 \\ +2 \\ +3\end{array}\right]$ + +$Ax=b$ + +### Vectors + +column vector x (length of 3): + +$\left[\begin{array}{c} +x_{1} \\ +x_{2} \\ +x_{3}\end{array}\right]$ + +row vector y (length of 3): + +$\left[\begin{array}{ccc} +y_{1}~ +y_{2}~ +y_{3}\end{array}\right]$ + +vector of length N: + +$\left[\begin{array}{c} +x_{1} \\ +x_{2} \\ +\vdots \\ +x_{N}\end{array}\right]$ + +The $i^{th}$ element of x is $x_{i}$ + + +```octave +x=[1:10] +``` + + x = + + 1 2 3 4 5 6 7 8 9 10 + + + + +```octave +x' +``` + + ans = + + 1 + 2 + 3 + 4 + 5 + 6 + 7 + 8 + 9 + 10 + + + +### Matrices + +Matrix A is 3x3: + +$A=\left[ \begin{array}{ccc} +5 & 3 & 0 \\ +1 & 2 & 3 \\ +1 & 1 & 1 \end{array} \right]$ + +elements in the matrix are denoted $A_{row~column}$, $A_{23}=3$ + +In general, matrix, B, can be any size, $M \times N$ (M-rows and N-columns): + +$B=\left[ \begin{array}{cccc} +B_{11} & B_{12} &...& B_{1N} \\ +B_{21} & B_{22} &...& B_{2N} \\ +\vdots & \vdots &\ddots& \vdots \\ +B_{M1} & B_{M2} &...& B_{MN}\end{array} \right]$ + + +```octave +A=[5,3,0;1,2,3;1,1,1] +``` + + A = + + 5 3 0 + 1 2 3 + 1 1 1 + + + +## Multiplication + +A column vector is a $1\times N$ matrix and a row vector is a $M\times 1$ matrix + +**Multiplying matrices is not commutative** + +$A B \neq B A$ + +Inner dimensions must agree, output is outer dimensions. + +A is $M1 \times N1$ and B is $M2 \times N2$ + +C=AB + +Therefore N1=M2 and C is $M1 \times N2$ + +If $C'=BA$, then N2=M1 and C is $M2 \times N1$ + + +e.g. +$A=\left[ \begin{array}{cc} +5 & 3 & 0 \\ +1 & 2 & 3 \end{array} \right]$ + +$B=\left[ \begin{array}{cccc} +1 & 2 & 3 & 4 \\ +5 & 6 & 7 & 8 \\ +9 & 10 & 11 & 12 \end{array} \right]$ + +C=AB + +$[2\times 4] = [2 \times 3][3 \times 4]$ + +The rule for multiplying matrices, A and B, is + +$C_{ij} = \sum_{k=1}^{n}A_{ik}B_{kj}$ + +In the previous example, + +$C_{11} = A_{11}B_{11}+A_{12}B_{21}+A_{13}B_{31} = 5*1+3*5+0*9=20$ + + +Multiplication is associative: + +$(AB)C = A(BC)$ + +and distributive: + +$A(B+C)=AB+AC$ + + +```octave +A=[5,3,0;1,2,3] +B=[1,2,3,4;5,6,7,8;9,10,11,12] +``` + + A = + + 5 3 0 + 1 2 3 + + B = + + 1 2 3 4 + 5 6 7 8 + 9 10 11 12 + + + + +```octave +C=zeros(2,4) +C(1,1)=A(1,1)*B(1,1)+A(1,2)*B(2,1)+A(1,3)*B(3,1); +C(1,2)=A(1,1)*B(1,2)+A(1,2)*B(2,2)+A(1,3)*B(3,2); + +C=A*B +``` + + C = + + 0 0 0 0 + 0 0 0 0 + + C = + + 20 28 36 44 + 38 44 50 56 + + + + +```octave +Cp=B*A +``` + + error: operator *: nonconformant arguments (op1 is 3x4, op2 is 2x3) + + +## Representation of a problem with Matrices and Vectors + +If you have a set of known output, $y_{1},~y_{2},~...y_{N}$ and a set of equations that +relate unknown inputs, $x_{1},~x_{2},~...x_{N}$, then these can be written in a vector +matrix format as: + +$y=Ax=\left[\begin{array}{c} +y_{1} \\ +y_{2} \\ +\vdots \\ +y_{N}\end{array}\right] += +A\left[\begin{array}{c} +x_{1} \\ +x_{2} \\ +\vdots \\ +x_{N}\end{array}\right]$ + +$A= +\left[\begin{array}{cccc} +| & | & & | \\ +a_{1} & a_{2} & ... & a_{N} \\ +| & | & & | \end{array}\right]$ + +or + +$y = a_{1}x_{1} + a_{2}x_{2} +...+a_{N}x_{N}$ + +where each $a_{i}$ is a column vector and $x_{i}$ is a scalar taken from the $i^{th}$ +element of x. + +Consider the following problem, 4 masses are connected in series to 4 springs with K=10 N/m. What are the final positions of the masses? + +![Springs-masses](mass_springs.svg) + +The masses haves the following amounts, 1, 2, 3, and 4 kg for masses 1-4. Using a FBD for each mass: + +$m_{1}g+k(x_{2}-x_{1})-kx_{1}=0$ + +$m_{2}g+k(x_{3}-x_{2})-k(x_{2}-x_{1})=0$ + +$m_{3}g+k(x_{4}-x_{3})-k(x_{3}-x_{2})=0$ + +$m_{4}g-k(x_{4}-x_{3})=0$ + +in matrix form: + +$\left[ \begin{array}{cccc} +2k & -k & 0 & 0 \\ +-k & 2k & -k & 0 \\ +0 & -k & 2k & -k \\ +0 & 0 & -k & k \end{array} \right] +\left[ \begin{array}{c} +x_{1} \\ +x_{2} \\ +x_{3} \\ +x_{4} \end{array} \right]= +\left[ \begin{array}{c} +m_{1}g \\ +m_{2}g \\ +m_{3}g \\ +m_{4}g \end{array} \right]$ + + +```octave +k=10; % N/m +m1=1; % kg +m2=2; +m3=3; +m4=4; +g=9.81; % m/s^2 +K=[2*k -k 0 0; -k 2*k -k 0; 0 -k 2*k -k; 0 0 -k k] +y=[m1*g;m2*g;m3*g;m4*g] + +x=K\y +``` + + K = + + 20 -10 0 0 + -10 20 -10 0 + 0 -10 20 -10 + 0 0 -10 10 + + y = + + 9.8100 + 19.6200 + 29.4300 + 39.2400 + + x = + + 9.8100 + 18.6390 + 25.5060 + 29.4300 + + + +## Matrix Operations + +Identity matrix `eye(M,N)` **Assume M=N unless specfied** + +$I_{ij} =1$ if $i=j$ + +$I_{ij} =0$ if $i\neq j$ + +AI=A=IA + +Diagonal matrix, D, is similar to the identity matrix, but each diagonal element can be any +number + +$D_{ij} =d_{i}$ if $i=j$ + +$D_{ij} =0$ if $i\neq j$ + +### Transpose + +The transpose of a matrix changes the rows -> columns and columns-> rows + +$(A^{T})_{ij} = A_{ji}$ + +for diagonal matrices, $D^{T}=D$ + +in general, the transpose has the following qualities: + +1. $(A^{T})^{T}=A$ + +2. $(AB)^{T} = B^{T}A^{T}$ + +3. $(A+B)^{T} = A^{T} +B^{T}$ + +All matrices have a symmetric and antisymmetric part: + +$A= 1/2(A+A^{T}) +1/2(A-A^{T})$ + +If $A=A^{T}$, then A is symmetric, if $A=-A^{T}$, then A is antisymmetric. + + +```octave +(A*B)' + +B'*A' % if this wasnt true, then inner dimensions wouldn;t match + +A'*B' == (A*B)' +``` + + ans = + + 20 38 + 28 44 + 36 50 + 44 56 + + ans = + + 20 38 + 28 44 + 36 50 + 44 56 + + error: operator *: nonconformant arguments (op1 is 3x2, op2 is 4x3) + + +## Square matrix operations + +### Trace + +The trace of a square matrix is the sum of the diagonal elements + +$tr~A=\sum_{i=1}^{N}A_{ii}$ + +The trace is invariant, meaning that if you change the basis through a rotation, then the +trace remains the same. + +It also has the following properties: + +1. $tr~A=tr~A^{T}$ + +2. $tr~A+tr~B=tr(A+B)$ + +3. $tr(kA)=k tr~A$ + +4. $tr(AB)=tr(BA)$ + + +```octave +id_m=eye(3) +trace(id_m) +``` + + id_m = + + Diagonal Matrix + + 1 0 0 + 0 1 0 + 0 0 1 + + ans = 3 + + +### Inverse + +The inverse of a square matrix, $A^{-1}$ is defined such that + +$A^{-1}A=I=AA^{-1}$ + +Not all square matrices have an inverse, they can be *singular* or *non-invertible* + +The inverse has the following properties: + +1. $(A^{-1})^{-1}=A$ + +2. $(AB)^{-1}=B^{-1}A^{-1}$ + +3. $(A^{-1})^{T}=(A^{T})^{-1}$ + + +```octave +A=rand(3,3) +``` + + A = + + 0.5762106 0.3533174 0.7172134 + 0.7061664 0.4863733 0.9423436 + 0.4255961 0.0016613 0.3561407 + + + + +```octave +Ainv=inv(A) +``` + + Ainv = + + 41.5613 -30.1783 -3.8467 + 36.2130 -24.2201 -8.8415 + -49.8356 36.1767 7.4460 + + + + +```octave +B=rand(3,3) +``` + + B = + + 0.524529 0.470856 0.708116 + 0.084491 0.730986 0.528292 + 0.303545 0.782522 0.389282 + + + + +```octave +inv(A*B) +inv(B)*inv(A) + +inv(A') + +inv(A)' +``` + + ans = + + -182.185 125.738 40.598 + -133.512 97.116 17.079 + 282.422 -200.333 -46.861 + + ans = + + -182.185 125.738 40.598 + -133.512 97.116 17.079 + 282.422 -200.333 -46.861 + + ans = + + 41.5613 36.2130 -49.8356 + -30.1783 -24.2201 36.1767 + -3.8467 -8.8415 7.4460 + + ans = + + 41.5613 36.2130 -49.8356 + -30.1783 -24.2201 36.1767 + -3.8467 -8.8415 7.4460 + + + +### Orthogonal Matrices + +Vectors are *orthogonal* if $x^{T}$ y=0, and a vector is *normalized* if $||x||_{2}=1$. A +square matrix is *orthogonal* if all its column vectors are orthogonal to each other and +normalized. The column vectors are then called *orthonormal* and the following results + +$U^{T}U=I=UU^{T}$ + +and + +$||Ux||_{2}=||x||_{2}$ + +### Determinant + +The **determinant** of a matrix has 3 properties + +1. The determinant of the identity matrix is one, $|I|=1$ + +2. If you multiply a single row by scalar $t$ then the determinant is $t|A|$: + +$t|A|=\left[ \begin{array}{cccc} +tA_{11} & tA_{12} &...& tA_{1N} \\ +A_{21} & A_{22} &...& A_{2N} \\ +\vdots & \vdots &\ddots& \vdots \\ +A_{M1} & A_{M2} &...& A_{MN}\end{array} \right]$ + +3. If you switch 2 rows, the determinant changes sign: + + +$-|A|=\left[ \begin{array}{cccc} +A_{21} & A_{22} &...& A_{2N} \\ +A_{11} & A_{12} &...& A_{1N} \\ +\vdots & \vdots &\ddots& \vdots \\ +A_{M1} & A_{M2} &...& A_{MN}\end{array} \right]$ + +4. inverse of the determinant is the determinant of the inverse: + +$|A^{-1}|=\frac{1}{|A|}=|A|^{-1}$ + +For a $2\times2$ matrix, + +$|A|=\left|\left[ \begin{array}{cc} +A_{11} & A_{12} \\ +A_{21} & A_{22} \\ +\end{array} \right]\right| = A_{11}A_{22}-A_{21}A_{12}$ + +For a $3\times3$ matrix, + +$|A|=\left|\left[ \begin{array}{ccc} +A_{11} & A_{12} & A_{13} \\ +A_{21} & A_{22} & A_{23} \\ +A_{31} & A_{32} & A_{33}\end{array} \right]\right|=$ + +$A_{11}A_{22}A_{33}+A_{12}A_{23}A_{31}+A_{13}A_{21}A_{32} +-A_{31}A_{22}A_{13}-A_{32}A_{23}A_{11}-A_{33}A_{21}A_{12}$ + +For larger matrices, the determinant is + +$|A|= +\frac{1}{n!} +\sum_{i_{r}}^{N} +\sum_{j_{r}}^{N} +\epsilon_{i_{1}...i_{N}} +\epsilon_{j_{1}...j_{N}} +A_{i_{1}j_{1}} +A_{i_{2}j_{2}} +... +A_{i_{N}j_{N}}$ + +Where the Levi-Cevita symbol is $\epsilon_{i_{1}i_{2}...i_{N}} = 1$ if it is an even +permutation and $\epsilon_{i_{1}i_{2}...i_{N}} = -1$ if it is an odd permutation and +$\epsilon_{i_{1}i_{2}...i_{N}} = 0$ if $i_{n}=i_{m}$ e.g. $i_{1}=i_{7}$ + +![Levi-Cevita rule for even and odd permutations for an $N\times N$ +determinant](even_odd_perm.svg) + +So a $4\times4$ matrix determinant, + +$|A|=\left|\left[ \begin{array}{cccc} +A_{11} & A_{12} & A_{13} & A_{14} \\ +A_{21} & A_{22} & A_{23} & A_{24} \\ +A_{31} & A_{32} & A_{33} & A_{34} \\ +A_{41} & A_{42} & A_{43} & A_{44} \end{array} \right]\right|$ + +$=\epsilon_{1234}\epsilon_{1234}A_{11}A_{22}A_{33}A_{44}+ +\epsilon_{2341}\epsilon_{1234}A_{21}A_{32}A_{43}A_{14}+ +\epsilon_{3412}\epsilon_{1234}A_{31}A_{42}A_{13}A_{24}+ +\epsilon_{4123}\epsilon_{1234}A_{41}A_{12}A_{23}A_{34}+... +\epsilon_{4321}\epsilon_{4321}A_{44}A_{33}A_{22}A_{11}$ + +### Special Case for determinants + +The determinant of a diagonal matrix $|D|=D_{11}D_{22}D_{33}...D_{NN}$. + +Similarly, if a matrix is upper triangular (so all values of $A_{ij}=0$ when $j + + + + + + + + + image/svg+xml + + + + + + + + + + + + + + m1 + m2 + m3 + + + m4 + + + + + + + + +