3/1/17

HW4/README.md

@@ -0,0 +1,37 @@
+# Homework #3
+## due 3/1/17 by 11:59pm
+
+
+1. Use your repository 'roots_and_optimization'. Document all the HW4 work under the
+heading `# Homework #4` in your `README.md` file
+
+    a. Create a function called 'collar_potential_energy' that computes the total
+    potential energy of a collar connected to a spring and sliding on a rod. As shown in
+    the figure given a position, xc, and angle, theta:
+
+    ![Collar-mass on an inclined rod](collar_mass.png)
+
+    The spring is unstretched when x_C=0.5 m. The potential energy due to gravity is:
+
+    PE_g=m x_C\*g\*sin(theta)
+
+    where m=0.5 kg, and g is the acceleration due to gravity, 
+
+    and the potential energy due to the spring is:
+
+    PE_s=1/2\*K \*(DL)^2
+
+    where DL = 0.5 - sqrt(0.5^2+(0.5-x_C)^2) and K=30 N/m.
+
+    b. Use the `goldmin.m` function to solve for the minimum potential energy at xc when
+    theta=0. *create an anonymous function with `@(x) collar_potential_energy(x,theta)` in
+    the input for goldmin. Be sure to include the script that solves for xc*
+
+    c. Create a for-loop that solves for the minimum potential energy position, xc, at a
+    given angle, theta, for theta = 0..90 degrees. 
+
+    d. Include a plot of xc vs theta. `plot(theta,xc)` with
+
+    `![Steady-state position of collar on rod at angle theta](plot.png)`
+
+3. Commit your changes to your repository. Sync your local repository with github.

README.md

@@ -79,10 +79,10 @@ general, I will not post homework solutions.
 |   |2/16|8|Linear Algebra|
 |6|2/21|9|Linear systems: Gauss elimination|
 |   |2/23|10|Linear Systems: LU factorization|
-|7|2/28||Midterm Review|
-|   |3/2||Midterm|
-|8|3/7|11|Linear Systems: Error analysis|
-|   |3/9|13|Eigenvalues|
+|7|2/28|11|Linear Systems: Error analysis|
+|   |3/2|12|Eigenvalues|
+|8|3/7|1-10 |Midterm Review|
+|   |3/9|1-10|Midterm|
 |9|3/14| N/A |Spring Break!|
 |   |3/16| N/A |Spring Break!|
 |10|3/21|12|Linear Systems: Iterative methods|

lecture_10/.ipynb_checkpoints/lecture_10-checkpoint.ipynb

@@ -0,0 +1,6 @@
+{
+ "cells": [],
+ "metadata": {},
+ "nbformat": 4,
+ "nbformat_minor": 2
+}

lecture_10/GaussNaive.m

@@ -0,0 +1,25 @@
+function [x,Aug] = GaussNaive(A,y)
+% 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+1;
+Aug = [A y];
+% forward elimination
+for k = 1:n-1
+  for i = k+1:n
+    factor = Aug(i,k)/Aug(k,k);
+    Aug(i,k:nb) = Aug(i,k:nb)-factor*Aug(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

lecture_10/GaussPivot.m

@@ -0,0 +1,33 @@
+function [x,Aug,npivots] = GaussPivot(A,b)
+% GaussPivot: Gauss elimination pivoting
+%   x = GaussPivot(A,b): Gauss elimination with pivoting.
+% input:
+%   A = coefficient matrix
+%   b = right hand side vector
+% output:
+%   x = solution vector
+[m,n]=size(A);
+if m~=n, error('Matrix A must be square'); end
+nb=n+1;
+Aug=[A b];
+npivots=0; % initially no pivots used
+% forward elimination
+for k = 1:n-1
+  % partial pivoting
+  [big,i]=max(abs(Aug(k:n,k)));
+  ipr=i+k-1;
+  if ipr~=k
+    npivots=npivots+1; % if the max is not the current index ipr, pivot count
+    Aug([k,ipr],:)=Aug([ipr,k],:);
+  end
+  for i = k+1:n
+    factor=Aug(i,k)/Aug(k,k);
+    Aug(i,k:nb)=Aug(i,k:nb)-factor*Aug(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