Lab 2 - Matrix multiplication, random numbers, plotting, program flow

Contents

1 Revisiting inner products and matrix multiplication

%Inner products
rowVector = [1 2 3 4];
columnVector = [5; 6; 7; 8];
innerProduct = rowVector*columnVector;

%Matrix multiplication
A = [1 2;3 4];
B = [5 6;7 8];
C = A*B
D = B*A
%Convince yourself that matrix multiplication is in fact *not* commutative,
%unless in very special cases

%Element-wise multiplication has be specified by a dot. MATLAB - for
%historical reasons assumes you want to do linear algebra by default.
E = A.*B
E = B.*A
%Good news: Your elementary school education is still valid. "Normal"
%plain vanilla element-wise multiplication commutes!

C =

    19    22
    43    50


D =

    23    34
    31    46


E =

     5    12
    21    32


E =

     5    12
    21    32

2 Introducing loops (repeats) by redoing inner product and matrix multiplication yourself

%1 Measure the number of elements of a vector. As they must the same number
%in both row and column vector, it is ok to measure just one. Later, when
%writing error-checking code, we can check whether they match before doing
%anything.
%2 Preallocate memory to speed up assignment of products to resulting
%product vector. We will use this vector to capture the products we make.
%In general, MATLAB loses a *lot* of time if you don't preallocate.
%3 Go through all the elements of the vectors, one by one. What ii is
%changes every time you go through the loop.
%4 For each element, multiply the one from the row vector with the
%corresponding one from the column vector and assign it to the right place
%in the products vector, replacing the 0s we preallocated
%5 Once we done looping, sum up all the products
numElements = length(rowVector) %1
products = zeros(numElements,1); %2
for ii = 1:numElements %3
  products(ii) =  rowVector(ii).*columnVector(ii); %4
end
innerProduct = sum(products) %5

%Exercise: Try different row and column vectors to convince yourself that
%this works in general (compare the in-built function to your loop).

%Let's do the same with matrix multiplication!
%1 Measure the number of rows of matrix A by calling size over dimension 1
%2 Measure the number of columns of matrix B by calling size over
%dimension2
%3 Preallocate a matrix with the right number of rows and columns with 0s.
%In the future, we will preallocate with nans because if we have a 0 in the
%resulting matrix, we don't know if one element wasn't assigned
%(overwritten) or whether the 0 is genuine. We don't hardcore - here or in
%general because we don't want this to work just for these specific
%matrices.
%4 Go through all the rows - you can call the counter variables anything
%you want. Make use of this, you don't have to call them i or j, like
%mathematicians do. Make it something good.
%5 Go through all the columns
%6 Replace the suitable zero in matrix C with the correct inner product
%7 Show C at the end
numRows = size(A,1) %1
numColumns = size(B,2) %2
C = zeros(numRows,numColumns); %3
for rr = 1:numRows %4
    for cc = 1:numColumns %5
        C(rr,cc) = A(rr,:)*B(:,cc); %6
    end
end
C %7

%Exercise: Try this for 2 other matrices and compare outcome with built-in
%matrix multiplication operator

numElements =

     4


innerProduct =

    70


numRows =

     2


numColumns =

     2


C =

    19    22
    43    50

3 Random numbers and plotting thereof and introducing conditionals

%A comment on comment pointers: If the comment exceeds one line, use a
%comment pointer. Like here.

%1 Show graph. Some versions of MATLAB and operating systems have this on by default
%If it doesn't the graph, it will still create it, just not show it. shg
%brings the figure you just created to the front of the screen.

A = randn(10000,2); %Randomly draw 10,000 numbers from a normal distribution
histogram(A,50); %Draw a histogram with 50 bins
shg %1

randSwitch = 1 %We decide which distribution to draw from, depending on state of switch
numMany = 100000; %This is how many numbers we want

%If statements do everything in order from if to end, but once a true
%statement is encountered, it is executed and the conditional exited.
if randSwitch == 1 %Draw randomly from a normal distribution
    A = randn(numMany,2); %2d numbers
elseif randSwitch == 2 %Draw randomly from a uniform distribution (0 1)
    A = rand(numMany,2); %2d numbers
elseif randSwitch == 3 %Draw random integers from 1 to 100
    A = randi(100,[numMany, 1]); %2d integers as a stack numMany high
else %Don't do anything
end

histogram(A,199)
shg
%xlim([-5 5]) %Make axis limits explicit. In the case of the normal distribution, something like -4 to plus 4 or -5 to 5 is reasonable because the dropoff is exponential.
xlim([min(A(:)) max(A(:))])

randSwitch =

     1

4 Introducing while statements and mixing if statements into loops

%Let's find 10 orthogonal vectors by random trial and error
orthoCounter = 1 ; %This is our counter of orthogonal vectors
loopCounter = 0; %This counts the number of loops it took to get them
while orthoCounter < 11 % We don't know how many tries we need, so while is suitable
    A = randn(1,2); %Draw a 2d row vector randomly from a normal distribution
    B = randn(2,1); %Draw a 2d column vector randomly from a normal distribution
    testThem = A*B; %Take the inner product
    loopCounter = loopCounter + 1;
    if abs(testThem) < 0.000001 %Ideally, we would IPs that are exactly 0, due to MATLABs numerical precision, this will never happen, so we need to check for close enough to zero
        disp('Orthogonal!') %Expression of joy
        rowContainer(orthoCounter,:) = A; %Capture the row vector
        columnContainer (:,orthoCounter) = B; %Capture the corresponding column vector. We should initialize these
        orthoCounter = orthoCounter + 1; %Increment the index
    end %Close the if statement #1

end %Close the while loop
%The reason I'm showing the next command is to give an example of how to
%mix characters and numbers in the same statement.
disp(['It took ', num2str(loopCounter), ' repeats to find 10 orthogonal vectors'])

Orthogonal!
Orthogonal!
Orthogonal!
Orthogonal!
Orthogonal!
Orthogonal!
Orthogonal!
Orthogonal!
Orthogonal!
Orthogonal!
It took 9552147 repeats to find 10 orthogonal vectors

5 Now that we found them, let's plot them

figure %Unless you call this explicitly and there is not one open, it will open one for you. Here, it is good to open one explictly
%If you have more than one figure you draw to, it helps to open new figures
for ii = 1:10
    plot([0 rowContainer(ii,1)],[0 rowContainer(ii,2)],'color','b') %Plot x-piece of row vector from origin, then y-piece in blue
    hold on %Put hold on, otherwise the 2nd vector will overwrite the first one
    plot([0 columnContainer(1,ii)], [0 columnContainer(2,ii)],'color','r') %Plot x-piece of column vector from origin, in red
    hold off %If I don't put the hold off again, the next iteration of the loop will all draw them on top of each other
    %Exercise: Play around with hold to see what will happen
        axis square
        axis equal %If axis units are not equal, it will look distorted
        pause %Only move on to the next loop on button press

end

6 Putting it all together - matrix multiplication, loops, random numbers

%We want to create a constellation of dots that represents some data or
%system
%We want to rotate it with a orthogonal matrix of unit vectors. These
%rotate the vectors without changing their lengths
close all %This closes all pre-existing figures
figure %Open a new figure
A = randn(1000,2); %Drawing 1000 random vectors in 2 dimensions from a normal distribution
theta = -1 / 180 * pi; %1 degree in radians
B = [cos(theta) sin(theta); -sin(theta) cos(theta)]; %One possible orthogonal unit vector
for ii = 1:359
    plot(A(:,1),A(:,2),'.','markersize',20); %Plot the endpoints of the vectors, not the vectors themselves. As dots
    xlim([-4 4]) %This will capture all dots from a normal distribution in all likelihood
    ylim([-4 4])
    pause(0.01) %With an argument, pause tells MATLAB how long to wait before moving on, in seconds

    A = A*B; %Updating A each iteration of the loop with the rotation from the ortogonal matrix
    shg
end

%This works beautifully. Two exercises:
%1 Make it rotate faster!
%2 Don't make it rotate, make it expand or contract.

Published with MATLAB® R2015a