Lab 9: Choosing parameters for physical simulations

Last lab we simulated a set of springs attached to two anchor points, and let the system evolve to a steady state. When that steady state is reached, the points should lie along a curve called a "catenary".

This curve has the form:

y = (a/2).*(exp(x/a) + exp(-x/a)
where exp is matlab for "e to the power of".
Task 1 make a function "catenary" that takes a parameter a, and a collection of points x and computes the y-value for those points x. Using commands like:
x = -3:0.1:3;
plot(x,catenary(0.5,x),'k'); hold on;
plot(x,catenary(1,x),'b');
plot(x,catenary(2,x),'r');  hold off;
Replicate the figure shown on the wiki page for catenary, and lookup the help for the "axis" command to make the plot have the same axes.
Task 2 When you finished lab 8, you had 2 matrices x,y that had the x,y locations of a collections of stones held together by springs. These are in the shape of a catenary, but do not exactly align because they are not centered around the x axis, and they have a different y offset. Take those points and center them with the command
xp = x-0.5
The goal of the remainder of this lab is to introduce matlabs fmins function which will help us search for the yoffset and the parameter a that minimizes a function. To do this we need to provide a function that computes the error for a given set of parameters (and we need the parameters that we want matlab to search over to be first). So create:
function err = catError(a,yoffset,xp,y)
  yp = yoffset + catenary(a,xp);     % compute the estimated y values
  err = ?????    % compute the sum of squared error between the
                   estimated y values and the values given in variable y
now, by hand try different values of a and yoffset until you get an error of less than 10. It may help you to write a helper routine to plot the points xp, y and the estimate catenary to see what the values are. (hints, y should be a small negative number and a should b between 0 and 1).

Task 3 We're almost ready to plug this into matlab, but the fmins function requires that the error function takes a vector to be optimized over so re-write your "catError" function use a 2 element vector P to replace a and yoffset.

function err = catError(P,xp,y)
  yp = P(2) + catenary(P(1),xp);     % compute the estimated y values
  err = ?????    % compute the sum of squared error between the
                   estimated y values and the values given in variable y
then you can try the following lines of code:
options = optimset('Display','iter');
pout = fminsearch('catError',[0.1 -1],options,xp,y);
Which follows the "simplex" method to minimize a function, and prints out a bunch of information including the value of the function you're trying to minimize. To see how good you did, you can display these results by using "pout" to reconstruct the best catenary, and then plotting the comparison.
cc = pout(2) + catenary(pout(1),xp);
plot(xp,y,'o-'); hold on;
plot(xp,cc,'r+-');
Show this plot the TA's to be checked off