mrw13002 · mrw13002 · Mar 7, 2017 · Mar 2, 2017 · Mar 2, 2017
diff --git a/lecture_12/lecture_12.ipynb b/lecture_12/lecture_12.ipynb
@@ -2,7 +2,7 @@
  "cells": [
   {
    "cell_type": "code",
-   "execution_count": 27,
+   "execution_count": 2,
    "metadata": {
     "collapsed": true
    },
@@ -13,7 +13,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 28,
+   "execution_count": 3,
    "metadata": {
     "collapsed": true
    },
@@ -702,7 +702,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 71,
+   "execution_count": 4,
    "metadata": {
     "collapsed": false
    },
@@ -828,7 +828,7 @@
        "\t</g>\n",
        "</g>\n",
        "<g color=\"black\" fill=\"none\" stroke=\"currentColor\" stroke-linecap=\"butt\" stroke-linejoin=\"miter\" stroke-width=\"3.00\">\n",
-       "\t<path d=\"M90.0,49.7 L143.8,49.7 M75.1,359.6 L79.8,366.2 L84.4,364.9 L89.1,365.8 L93.7,366.6 L98.4,366.2   L103.0,365.3 L107.7,365.3 L112.3,363.8 L117.0,363.2 L121.6,360.0 L126.2,347.8 L130.9,358.3 L135.5,362.7   L140.2,362.4 L144.8,361.4 L149.5,360.1 L154.1,323.0 L158.8,341.7 L163.4,341.7 L168.0,332.4 L172.7,338.1   L177.3,335.4 L182.0,342.0 L186.6,364.3 L191.3,371.9 L195.9,370.8 L200.6,371.3 L205.2,371.2 L209.9,370.2   L214.5,368.8 L219.1,369.7 L223.8,368.4 L228.4,364.2 L233.1,366.7 L237.7,366.9 L242.4,364.7 L247.0,365.3   L251.7,364.0 L256.3,364.6 L260.9,363.4 L265.6,362.2 L270.2,361.6 L274.9,361.6 L279.5,360.3 L284.2,358.3   L288.8,358.7 L293.5,359.0 L298.1,357.3 L302.8,354.8 L307.4,355.2 L312.0,355.4 L316.7,353.9 L321.3,353.1   L326.0,351.9 L330.6,352.8 L335.3,349.5 L339.9,350.0 L344.6,348.2 L349.2,349.1 L353.8,345.0 L358.5,346.2   L363.1,344.9 L367.8,346.9 L372.4,337.9 L377.1,341.0 L381.7,342.9 L386.4,341.2 L391.0,339.7 L395.7,339.0   L400.3,338.1 L404.9,337.3 L409.6,336.8 L414.2,336.3 L418.9,334.4 L423.5,330.6 L428.2,330.4 L432.8,332.9   L437.5,329.8 L442.1,330.6 L446.7,327.8 L451.4,326.2 L456.0,323.9 L460.7,323.0 L465.3,322.7 L470.0,321.2   L474.6,320.8 L479.3,312.8 L483.9,317.3 L488.6,317.1 L493.2,317.3 L497.8,315.7 L502.5,316.6 L507.1,315.7   L511.8,307.6 L516.4,315.5 L521.1,314.6 L525.7,311.8 L530.4,308.5 L535.0,310.2  \" stroke=\"rgb(  0,   0, 255)\"/></g>\n",
+       "\t<path d=\"M90.0,49.7 L143.8,49.7 M75.1,365.8 L79.8,371.5 L84.4,373.0 L89.1,373.5 L93.7,373.1 L98.4,372.6   L103.0,372.2 L107.7,371.9 L112.3,371.1 L117.0,371.3 L121.6,370.2 L126.2,370.6 L130.9,367.8 L135.5,368.4   L140.2,367.8 L144.8,366.9 L149.5,366.4 L154.1,332.1 L158.8,359.6 L163.4,361.2 L168.0,360.5 L172.7,360.3   L177.3,357.2 L182.0,357.6 L186.6,351.1 L191.3,355.9 L195.9,351.5 L200.6,353.7 L205.2,351.7 L209.9,351.3   L212.7,16.7 M216.4,16.7 L219.1,330.9 L223.8,327.4 L228.4,324.3 L233.1,326.2 L237.7,323.8 L242.4,322.5   L247.0,319.7 L251.7,316.8 L256.3,317.0 L260.9,312.9 L265.6,311.3 L270.2,304.0 L274.9,308.9 L279.5,306.3   L284.2,329.1 L288.8,327.3 L293.5,327.4 L298.1,323.6 L302.8,336.6 L307.4,336.1 L312.0,336.6 L316.7,333.9   L321.3,331.9 L326.0,330.4 L330.6,331.5 L335.3,326.5 L339.9,337.7 L344.6,335.9 L349.2,336.6 L353.8,333.3   L358.5,327.6 L363.1,331.2 L367.8,328.9 L372.4,327.4 L377.1,335.1 L381.7,333.5 L386.4,336.3 L391.0,332.8   L395.7,332.2 L400.3,331.1 L404.9,331.3 L409.6,328.0 L414.2,336.4 L418.9,331.6 L423.5,334.8 L428.2,333.9   L432.8,329.7 L437.5,329.6 L442.1,331.3 L446.7,325.6 L451.4,332.9 L456.0,332.4 L460.7,332.4 L465.3,328.0   L470.0,329.8 L474.6,322.5 L479.3,326.9 L483.9,325.8 L488.6,322.3 L493.2,322.5 L497.8,323.4 L502.5,322.1   L507.1,321.9 L511.8,319.4 L516.4,322.7 L521.1,318.5 L525.7,315.9 L530.4,315.0 L535.0,312.2  \" stroke=\"rgb(  0,   0, 255)\"/></g>\n",
        "\t</g>\n",
        "\t<g id=\"gnuplot_plot_2a\"><title>backslash</title>\n",
        "<g color=\"black\" fill=\"none\" stroke=\"currentColor\" stroke-linecap=\"butt\" stroke-linejoin=\"miter\" stroke-width=\"3.00\">\n",
@@ -837,7 +837,7 @@
        "\t</g>\n",
        "</g>\n",
        "<g color=\"black\" fill=\"none\" stroke=\"currentColor\" stroke-linecap=\"butt\" stroke-linejoin=\"miter\" stroke-width=\"3.00\">\n",
-       "\t<path d=\"M90.0,97.7 L143.8,97.7 M75.1,379.4 L79.8,372.8 L84.4,373.9 L89.1,368.0 L93.7,373.0 L98.4,373.0   L103.0,372.1 L107.7,372.1 L112.3,371.2 L117.0,370.6 L121.6,367.8 L126.2,366.0 L130.9,366.7 L135.5,369.9   L140.2,369.9 L144.8,369.3 L149.5,368.5 L154.1,354.8 L158.8,355.7 L163.4,355.4 L168.0,350.5 L172.7,354.1   L177.3,349.7 L182.0,334.4 L186.6,376.5 L191.3,376.1 L195.9,375.5 L200.6,375.8 L205.2,375.9 L209.9,375.0   L214.5,374.8 L219.1,374.6 L223.8,374.1 L228.4,370.8 L233.1,373.0 L237.7,373.0 L242.4,371.3 L247.0,371.9   L251.7,371.5 L256.3,371.7 L260.9,371.5 L265.6,370.6 L270.2,370.2 L274.9,369.9 L279.5,369.5 L284.2,369.1   L288.8,368.8 L293.5,368.7 L298.1,368.2 L302.8,367.5 L307.4,367.1 L312.0,366.9 L316.7,366.4 L321.3,366.0   L326.0,365.6 L330.6,365.8 L335.3,364.9 L339.9,364.5 L344.6,363.8 L349.2,363.8 L353.8,361.0 L358.5,362.9   L363.1,360.5 L367.8,363.2 L372.4,361.9 L377.1,357.7 L381.7,360.0 L386.4,359.6 L391.0,358.9 L395.7,358.6   L400.3,357.7 L404.9,356.1 L409.6,357.4 L414.2,356.8 L418.9,356.3 L423.5,357.7 L428.2,353.1 L432.8,355.2   L437.5,354.2 L442.1,352.8 L446.7,351.7 L451.4,350.8 L456.0,350.2 L460.7,350.6 L465.3,348.7 L470.0,348.7   L474.6,348.0 L479.3,343.4 L483.9,347.4 L488.6,345.2 L493.2,347.3 L497.8,347.8 L502.5,346.7 L507.1,346.2   L511.8,338.6 L516.4,347.8 L521.1,346.4 L525.7,344.3 L530.4,345.6 L535.0,347.1  \" stroke=\"rgb(  0, 128,   0)\"/></g>\n",
+       "\t<path d=\"M90.0,97.7 L143.8,97.7 M75.1,380.8 L79.8,375.5 L84.4,377.0 L89.1,377.2 L93.7,377.2 L98.4,376.9   L103.0,376.5 L107.7,376.5 L112.3,375.9 L117.0,376.1 L121.6,375.2 L126.2,374.8 L130.9,373.0 L135.5,373.7   L140.2,373.2 L144.8,372.4 L149.5,372.2 L154.1,368.0 L158.8,368.6 L163.4,368.4 L168.0,367.7 L172.7,367.6   L177.3,365.4 L182.0,366.4 L186.6,363.1 L191.3,364.9 L195.9,362.9 L200.6,363.2 L205.2,362.5 L209.9,362.2   L214.5,353.0 L219.1,351.7 L223.8,347.4 L228.4,347.3 L233.1,347.3 L237.7,345.8 L242.4,344.7 L247.0,344.1   L251.7,342.3 L256.3,342.1 L260.9,340.7 L265.6,339.2 L270.2,333.7 L274.9,336.8 L279.5,351.7 L284.2,350.9   L288.8,350.2 L292.9,16.7 M294.0,16.7 L298.1,357.4 L302.8,355.5 L307.4,355.7 L312.0,355.5 L316.7,353.7   L321.3,354.2 L326.0,352.4 L330.6,353.4 L335.3,357.9 L339.9,357.6 L344.6,356.6 L349.2,356.9 L353.8,355.9   L358.5,353.2 L363.1,355.0 L367.8,354.8 L372.4,358.5 L377.1,355.7 L381.7,353.0 L386.4,355.5 L391.0,354.6   L395.7,353.7 L400.3,353.4 L404.9,353.7 L409.6,356.8 L414.2,356.5 L418.9,355.7 L423.5,355.7 L428.2,355.0   L432.8,352.3 L437.5,353.7 L442.1,352.4 L446.7,355.3 L451.4,355.2 L456.0,354.1 L460.7,353.7 L465.3,350.0   L470.0,353.2 L474.6,352.0 L479.3,350.0 L483.9,351.7 L488.6,348.9 L493.2,345.8 L497.8,350.9 L502.5,349.7   L507.1,349.5 L511.8,348.9 L516.4,349.5 L521.1,348.8 L525.7,346.3 L530.4,346.7 L535.0,343.8  \" stroke=\"rgb(  0, 128,   0)\"/></g>\n",
        "\t</g>\n",
        "\t<g id=\"gnuplot_plot_3a\"><title>multiplication</title>\n",
        "<g color=\"black\" fill=\"none\" stroke=\"currentColor\" stroke-linecap=\"butt\" stroke-linejoin=\"miter\" stroke-width=\"3.00\">\n",
@@ -846,7 +846,7 @@
        "\t</g>\n",
        "</g>\n",
        "<g color=\"black\" fill=\"none\" stroke=\"currentColor\" stroke-linecap=\"butt\" stroke-linejoin=\"miter\" stroke-width=\"3.00\">\n",
-       "\t<path d=\"M90.0,145.7 L143.8,145.7 M75.1,380.5 L79.8,379.4 L84.4,379.6 L89.1,379.8 L93.7,379.6 L98.4,379.8   L103.0,379.8 L107.7,379.8 L112.3,379.8 L117.0,379.8 L121.6,378.1 L126.2,379.4 L130.9,380.1 L135.5,380.3   L140.2,380.3 L144.8,380.3 L149.5,380.3 L154.1,378.9 L158.8,378.8 L163.4,378.9 L168.0,378.7 L172.7,378.5   L177.3,379.6 L182.0,379.4 L186.6,382.7 L191.3,382.7 L195.9,382.7 L200.6,382.7 L205.2,382.7 L209.9,382.7   L214.5,382.9 L219.1,382.9 L223.8,382.9 L228.4,382.7 L233.1,382.7 L237.7,382.7 L242.4,382.7 L247.0,382.9   L251.7,382.7 L256.3,382.7 L260.9,382.6 L265.6,382.7 L270.2,382.7 L274.9,382.7 L279.5,382.7 L284.2,382.7   L288.8,382.5 L293.5,382.7 L298.1,382.7 L302.8,382.6 L307.4,382.5 L312.0,382.5 L316.7,382.5 L321.3,382.7   L326.0,382.7 L330.6,382.7 L335.3,382.6 L339.9,382.7 L344.6,382.7 L349.2,382.6 L353.8,382.7 L358.5,382.6   L363.1,382.7 L367.8,382.5 L372.4,382.5 L377.1,382.7 L381.7,382.5 L386.4,382.7 L391.0,382.5 L395.7,382.5   L400.3,382.6 L404.9,382.7 L409.6,382.6 L414.2,382.5 L418.9,382.6 L423.5,382.5 L428.2,382.5 L432.8,382.4   L437.5,382.5 L442.1,382.6 L446.7,382.6 L451.4,382.5 L456.0,382.4 L460.7,382.3 L465.3,382.5 L470.0,382.3   L474.6,382.6 L479.3,382.3 L483.9,382.4 L488.6,382.3 L493.2,382.3 L497.8,382.5 L502.5,382.3 L507.1,382.6   L511.8,382.4 L516.4,379.8 L521.1,382.3 L525.7,382.3 L530.4,382.5 L535.0,382.5  \" stroke=\"rgb(255,   0,   0)\"/></g>\n",
+       "\t<path d=\"M90.0,145.7 L143.8,145.7 M75.1,381.6 L79.8,381.1 L84.4,381.2 L89.1,381.2 L93.7,381.4 L98.4,381.2   L103.0,381.2 L107.7,381.4 L112.3,381.4 L117.0,381.2 L121.6,381.2 L126.2,381.4 L130.9,381.4 L135.5,381.2   L140.2,381.2 L144.8,381.2 L149.5,381.2 L154.1,380.9 L158.8,381.2 L163.4,381.2 L168.0,381.1 L172.7,381.2   L177.3,381.2 L182.0,381.2 L186.6,381.1 L191.3,381.1 L195.9,381.2 L200.6,378.1 L205.2,381.2 L209.9,381.1   L214.5,379.4 L219.1,380.0 L223.8,379.9 L228.4,380.1 L233.1,380.0 L237.7,380.1 L242.4,380.0 L247.0,379.9   L251.7,380.0 L256.3,379.8 L260.9,379.8 L265.6,380.0 L270.2,379.8 L274.9,379.8 L279.5,381.2 L284.2,381.2   L288.8,381.1 L293.5,381.1 L298.1,382.0 L302.8,381.8 L307.4,381.8 L312.0,381.8 L316.7,381.8 L321.3,381.8   L326.0,381.8 L330.6,381.8 L335.3,382.2 L339.9,382.2 L344.6,382.2 L349.2,382.2 L353.8,382.2 L358.5,381.9   L363.1,382.2 L367.8,382.2 L372.4,382.3 L377.1,382.5 L381.7,382.3 L386.4,382.4 L391.0,382.4 L395.7,382.3   L400.3,382.3 L404.9,381.6 L409.6,382.3 L414.2,382.5 L418.9,382.6 L423.5,382.6 L428.2,382.6 L432.8,382.5   L437.5,382.5 L442.1,382.5 L446.7,382.7 L451.4,382.7 L456.0,382.7 L460.7,382.5 L465.3,382.5 L470.0,382.6   L474.6,382.5 L479.3,382.7 L483.9,382.7 L488.6,382.6 L493.2,382.7 L497.8,382.6 L502.5,382.5 L507.1,382.6   L511.8,382.5 L516.4,382.4 L521.1,382.3 L525.7,382.5 L530.4,382.3 L535.0,382.3  \" stroke=\"rgb(255,   0,   0)\"/></g>\n",
        "\t</g>\n",
        "<g color=\"white\" fill=\"none\" stroke=\"rgb(255,   0,   0)\" stroke-linecap=\"butt\" stroke-linejoin=\"miter\" stroke-width=\"2.00\">\n",
        "</g>\n",

diff --git a/lecture_13/.GS_rel.m.swp b/lecture_13/.GS_rel.m.swp
diff --git a/lecture_13/.GaussSeidel.m.swp b/lecture_13/.GaussSeidel.m.swp
diff --git a/lecture_13/.Jacobi_rel.m.swp b/lecture_13/.Jacobi_rel.m.swp
diff --git a/...b_checkpoints/lecture_12-checkpoint.ipynb → ...b_checkpoints/lecture_13-checkpoint.ipynb b/...b_checkpoints/lecture_12-checkpoint.ipynb → ...b_checkpoints/lecture_13-checkpoint.ipynb
diff --git a/lecture_13/.lambda_fcn.m.swp b/lecture_13/.lambda_fcn.m.swp
diff --git a/lecture_13/GS_rel.m b/lecture_13/GS_rel.m
@@ -0,0 +1,36 @@
+function [x,ea,iter] = GS_rel(A,b,lambda,es,maxit)
+% GaussSeidel: Gauss Seidel method
+%   x = GaussSeidel(A,b): Gauss Seidel without relaxation
+% input:
+%   A = coefficient matrix
+%   b = right hand side vector
+%   es = stop criterion (default = 0.00001%)
+%   maxit = max iterations (default = 50)
+% output:
+%   x = solution vector
+if nargin<3,error('at least 2 input arguments required'),end
+if nargin<5|isempty(maxit),maxit=50;end
+if nargin<4|isempty(es),es=0.00001;end
+[m,n] = size(A);
+if m~=n, error('Matrix A must be square'); end
+C = A-diag(diag(A));
+x=zeros(n,1);
+for i = 1:n
+  C(i,1:n) = C(i,1:n)/A(i,i);
+end
+
+d = b./diag(A);
+
+iter = 0;
+while (1)
+  xold = x;
+  for i = 1:n
+    x(i) = d(i)-C(i,:)*x;
+    x(i) = lambda*x(i)+(1-lambda)*xold(i);
+    if x(i) ~= 0
+      ea(i) = abs((x(i) - xold(i))/x(i)) * 100;
+    end
+  end
+  iter = iter+1;
+  if max(ea)<=es | iter >= maxit, break, end
+end
diff --git a/lecture_13/GaussSeidel.m b/lecture_13/GaussSeidel.m
@@ -0,0 +1,35 @@
+function x = GaussSeidel(A,b,es,maxit)
+% GaussSeidel: Gauss Seidel method
+%   x = GaussSeidel(A,b): Gauss Seidel without relaxation
+% input:
+%   A = coefficient matrix
+%   b = right hand side vector
+%   es = stop criterion (default = 0.00001%)
+%   maxit = max iterations (default = 50)
+% output:
+%   x = solution vector
+if nargin<2,error('at least 2 input arguments required'),end
+if nargin<4|isempty(maxit),maxit=50;end
+if nargin<3|isempty(es),es=0.00001;end
+[m,n] = size(A);
+if m~=n, error('Matrix A must be square'); end
+C = A-diag(diag(A));
+x=zeros(n,1);
+for i = 1:n
+  C(i,1:n) = C(i,1:n)/A(i,i);
+end
+
+d = b./diag(A);
+
+iter = 0;
+while (1)
+  xold = x;
+  for i = 1:n
+    x(i) = d(i)-C(i,:)*x;
+    if x(i) ~= 0
+      ea(i) = abs((x(i) - xold(i))/x(i)) * 100;
+    end
+  end
+  iter = iter+1;
+  if max(ea)<=es | iter >= maxit, break, end
+end
diff --git a/lecture_13/Jacobi.m b/lecture_13/Jacobi.m
@@ -0,0 +1,39 @@
+function x = Jacobi(A,b,es,maxit)
+% GaussSeidel: Gauss Seidel method
+%   x = GaussSeidel(A,b): Gauss Seidel without relaxation
+% input:
+%   A = coefficient matrix
+%   b = right hand side vector
+%   es = stop criterion (default = 0.00001%)
+%   maxit = max iterations (default = 50)
+% output:
+%   x = solution vector
+if nargin<2,error('at least 2 input arguments required'),end
+if nargin<4|isempty(maxit),maxit=50;end
+if nargin<3|isempty(es),es=0.00001;end
+[m,n] = size(A);
+if m~=n, error('Matrix A must be square'); end
+C = A-diag(diag(A));
+x=zeros(n,1);
+for i = 1:n
+  C(i,1:n) = C(i,1:n)/A(i,i);
+end
+
+d = b./diag(A);
+
+iter = 0;
+while (1)
+  xold = x;
+  x = d-C*x;
+  % if any values of x are zero, we add 1 to denominator so error is well-behaved
+  i_zero=find(x==0);
+  i=find(x~=0);
+  if length(i_zero)>0
+    ea(i_zero)=abs((x-xold)./(1+x)*100);
+    ea(i) = abs((x(i) - xold(i))./x(i)) * 100;
+  else
+    ea = abs((x - xold)./x) * 100;
+  end
+  iter = iter+1;
+  if max(ea)<=es | iter >= maxit, break, end
+end
diff --git a/lecture_13/Jacobi_rel.m b/lecture_13/Jacobi_rel.m
@@ -0,0 +1,41 @@
+function [x,ea,iter]= Jacobi_rel(A,b,lambda,es,maxit)
+% GaussSeidel: Gauss Seidel method
+%   x = GaussSeidel(A,b): Gauss Seidel without relaxation
+% input:
+%   A = coefficient matrix
+%   b = right hand side vector
+%   es = stop criterion (default = 0.00001%)
+%   maxit = max iterations (default = 50)
+% output:
+%   x = solution vector
+if nargin<3,error('at least 2 input arguments required'),end
+if nargin<5|isempty(maxit),maxit=50;end
+if nargin<4|isempty(es),es=0.00001;end
+[m,n] = size(A);
+if m~=n, error('Matrix A must be square'); end
+C = A-diag(diag(A));
+x=zeros(n,1);
+for i = 1:n
+  C(i,1:n) = C(i,1:n)/A(i,i);
+end
+
+d = b./diag(A);
+
+iter = 0;
+while (1)
+  xold = x;
+  x = d-C*x;
+  % Add relaxation parameter lambda to current iteration
+  x = lambda*x+(1-lambda)*xold;
+  % if any values of x are zero, we add 1 to denominator so error is well-behaved
+  i_zero=find(x==0);
+  i=find(x~=0);
+  if length(i_zero)>0
+    ea(i_zero)=abs((x-xold)./(1+x)*100);
+    ea(i) = abs((x(i) - xold(i))./x(i)) * 100;
+  else
+    ea = abs((x - xold)./x) * 100;
+  end
+  iter = iter+1;
+  if max(ea)<=es | iter >= maxit, break, end
+end
diff --git a/lecture_13/LU_naive.m b/lecture_13/LU_naive.m
@@ -0,0 +1,27 @@
+function [L, U] = LU_naive(A)
+% GaussNaive: naive Gauss elimination
+%   x = GaussNaive(A,b): Gauss elimination without pivoting.
+% input:
+%   A = coefficient matrix
+%   y = right hand side vector
+% output:
+%   x = solution vector
+[m,n] = size(A);
+if m~=n, error('Matrix A must be square'); end
+nb = n;
+L=diag(ones(n,1));
+U=A;
+% forward elimination
+for k = 1:n-1
+  for i = k+1:n
+    fik = U(i,k)/U(k,k);
+    L(i,k)=fik;
+    U(i,k:nb) = U(i,k:nb)-fik*U(k,k:nb);
+  end
+end
+%% back substitution
+%x = zeros(n,1);
+%x(n) = Aug(n,nb)/Aug(n,n);
+%for i = n-1:-1:1
+%  x(i) = (Aug(i,nb)-Aug(i,i+1:n)*x(i+1:n))/Aug(i,i);
+%end
diff --git a/lecture_13/efficient_soln.png b/lecture_13/efficient_soln.png
diff --git a/lecture_13/gp_image_01.png b/lecture_13/gp_image_01.png
diff --git a/lecture_13/lambda_fcn.m b/lecture_13/lambda_fcn.m
@@ -0,0 +1,8 @@
+function iters = lambda_fcn(A,b,lambda)
+  % function to minimize the number of iterations for a given Ax=b solution
+  % using default Jacobi_rel parameters of es=0.00001 and maxit=50
+
+  [x,ea,iters]= Jacobi_rel(A,b,lambda,1e-8);
+end
+
+