Objectives

• Use Matlab to solve a system of equations.
• Write a function containing MATLAB code to implement part of MP1
• Use the 'diag' and 'length' and 'cumsum' Matlab functions.

References

• Lecture 4 on MATLAB  matrices.
• MATLAB Help: Using MATLAB  Programming and Data Types  M-File Programming
• Matlab by Pratap Chapters 3.1, 3.2, 5.1

 

Instructions

• In this lab, you will work individually. You can ask your group members for help but each person must submit their own answer sheet and write the find_xy.m function.
  
• Read the instructions for each part before filling out the answer sheet.

Part 1

Part 2: Programming the find_xy function.

In this lab we will be writing MATLAB code to implement the central part of the solution to the suspension bridge problem from MP1. Take time now to read through Part 1 and Part 2(shows the equation Dy =  wp) of MP1. The ultimate goal of this lab is to get MATLAB to solve a system of equations represented by a matrix equation similar to the one in MP1. To solve the system, we need to solve the matrix equation:
D*y = wp

where y is a nx1 column vector, D is an nxn matrix and wp is an nx1 column vector.

Your code should perform the calculations necessary to solve the matrix equation above to obtain the value of the vector y, which represents the values of the heights of the cables. The data characterizing the system are provided in the file mp1_data.m. You can download this file here. (Right click and save it to your home directory.)

The file contains the following variables that describe the known quantities of the system.

• n: The number of support cables.
• L: Length of the main cable.
• yL: the height of the left tower. (this is equal to y₀ in the set of equations below)
• yR: the height of the right tower. (this is equal to y_n+1 in the set of equation below)
• d: A row vector of length n+1 containing the horizontal distances between the left tower and the 1^st cable, between each pair of consecutive cables and between the right-most(n^th) cable and the right tower.
• W: A column vector of length n containing the weight on each cable.

      y_i-1*1/d_i + y_i*(-1/d_i + -1/d_i+1) + y_i+1 * 1/d_i+1 = W_i/ HF for i=1 .. n
 

These equations define a tri-diagonal system of linear equations for the heights y₁ to y_n corresponding to the weights W₁ to W_n.

For this lab, you will write a function containing MATLAB code called find_xy.m that you can execute from within MATLAB by typing "find_xy(HF)" at the MATLAB prompt. Of course you must assign a value for HF before you can type 'find_xy(HF)'.

1. Open Matlab and edit a new file called "find_xy.m" . Add the following lines to your find_xy.m file:

 

function[x,y] = find_xy(HF)
% Begin the file by writing a few lines of
% comments (lines that start with '%').
% Your comments should include your names, the date,
% and a brief  description of the program you are about to write.
mp1_data     % mp1_data.m is a script file.
% You can run a script file from inside a function.
% Running a script from inside the function
% find_xy make the values of the
% variables in mp1_data available to the function find_xy.

% Begin typing your  code after this comment.

 

For the next two steps you will go over the math behind the problem for MP1. Don't do any code for your Matlab function until step 3.



2. Using the equation for y_i given above, write out the equations for i = 1.. 5 (n = 5) on the answer sheet for question #1.
   
   HINTS: y₀ = yL (the left tower), y₆ = yR (the right tower). These aren't part of the vector y that you're calculating; vector y has only n elements, the n cables. yR and yL are given to you as data characterizing the system, so they can be treated as constants. All constants can be moved to the right side of the equation.

 
 
3. Using the equations from the previous step , answer question #2 by writing out the matrix equation D*y = wp (don't use Matlab notation) on the answer sheet. This will make it a lot easier to understand the rest of the lab.  See the MP1, at the end of Part 2 where it shows how you should write the matrix equation D*y= wp . The example shown in MP1 is for any general n. To answer question #2 use n = 5.



From this point on, you will be writing MATLAB code to work with vectors and matrices. Remember that all the code you write should go into your find_xy.m file. Although the mp1_data.m datafile sets n=5 and we wrote the equations for n=5 above, the Matlab code you now write should work for any n.

a) Now you should start putting together the pieces of the matrices in the equation above. We'll start with the D matrix on the left-hand side of the equation mentioned above. We will do this by creating three matrices, D_main, D_upper and D_lower so that
                D = D_main + D_upper + D_lower ;


  Don't type this in your Matlab function yet as we need to first create the three matrices  D_main, D_upper and D_lower.
 D_main is a diagonal matrix whose only non-zero values lie on the main diagonal and are   -(1/d_i + 1/d_i+1)  i = 1,...n  .

Answer the question #3 on the answer sheet.

D_upper is the upper diagonal matrix whose only non-zero values lie above the main diagonal and are  1/d_i  i = 2,...n  .
 
Answer the question #4 on the answer sheet.
 

D_lower is the lower diagonal matrix whose only non-zero values lie below the main diagonal and are  1/d_i  i = 2,...n  .

Answer the question #5 on the answer sheet.

Now we can type

    D = D_main + D_upper + D_lower ;
in you find_xy.m file.



b) We need to create the right-hand vector wp.  This vector is created using the weight vector W, the HF value and the values yL and yR, all of which are given to you. Except for the first and last elements of wp, all the others have the form:
 


wp_i = W_i / HF;

Answer question #6.

The first element is W₁ / HF - yL/d₁ .The last element is W_n / HF - yR/d_n+1
Therefore all we need to do is subtract yL/d₁  from the first element of wp and subtract yR/d_n+1 from the nth element of wp.
Answer question #7.



c) Since we know the left-hand side matrix D and the right-hand side vector wp, we can solve the system of equations D*y = wp.
Store the solution in a column vector called y. But since our function find_xy returns two row vectors x and y we need to transpose vector y.
Also y should contain the heights of the left tower yL and right tower yR, so we will need to amend y by adding these values.


Answer questions #8,#9 and #10.



d) Create a vector x to hold the corresponding x-coordinates of the towers and the cables, so that you can plot all of these points. The distances are measured from the left tower. They are derived from the d vector that you are given (in mp1_data.m). Your first value in this

x vector should be 0 (corresponding to the left tower), and the rest of them should cumulatively add up the distances between the left tower and the first cable, each pair of successive cables, and the last cable and the right tower (think cumsum).


Answer question #11. We’re finished with the function at this point, so you can save it.




Part 3 Testing the find_xy function

Now that you have solved the problem correctly, you should test your code.  After saving your code in 'find_xy.m' , then in the Matlab command window type:
>> HF = 15
>> [x,y] = find_xy(HF)
You should get the following values for x and y:

x =

  0  5.3931  10.3278  14.6056  20.0877  25.0048  30.1779
y =
  30.0000  23.4883  19.8001  18.1607  18.3307  20.4408  25.0000
If you don't get the above values, you must "debug" your program.( Note: we are displaying only 4 digits to the right of the decimal. Rounding may occur, so if you display more digits you answers may vary.)  If you get an error message, read it carefully. Some messages are cryptic and therefore useless however most contain useful information like the line number in the file find_xy where the error occurred. When this doesn't work ask your TA for help.

 

Congratulations, you are done!