MATLAB PROJECT 1

Please read the Instructions located on the Assignments page prior

to working on the Project.

BEGIN with creating Live Script Project1.

Note: All exercises in this project will be completed in the Live Script using the Live Editor.

Please refer to the MATLAB video that explains how to use the Live Script:

https://www.mathworks.com/videos/using-the-live-editor-117940.html?s_tid=srchtitle

The final script has to be generated by exporting the Live Script to PDF.

Each exercise has to begin with the line

Exercise#

You should also mark down the parts such as (a), (b), (c), and etc. This makes grading easier.

Important: we use the default format short for the numbers in all exercises unless it is

specified otherwise. We do not employ format rat since it may cause problems with

running the codes and displaying matrices in Live Script. If format long has been used,

please make sure to return to the default format in the next exercise.

.

EXERCISE1 (4 points) Difficulty: Easy

Jordan Block is a square n n × matrix ( n ≥ 2 ) with a scalar r on the main diagonal and 1’s on

the diagonal right above it. All other entries are zero.

**Create a function in a file that begins with

function J=jord(n,r)

which, when possible, produces an n n × Jordan matrix with a scalar r on its main diagonal.

**First, the function will verify whether n >1 and n is an integer number. You can use a builtin MATLAB function such as mod or fix .

**If n does not satisfies the above conditions, the code has to assign an empty matrix to J and

display a message that ‘Jordan Block cannot be built’- you can use disp, sprintf, or

fprintf commands.

**If n is an integer number and n >1, you will proceed with constructing a Jordan Block

matrix J as described at the beginning of this exercise.

After the function jord is created and saved in a file, we return to the Live Script.

**First, we set the format mode by typing

format

format compact

**Then, we print in the Live Script the created function jord by typing

type jord

**Next, we input

r=rand(1)

2

and run the function J=jord(n,r) for each of the following n to get either a non-empty

output J or an empty matrix and a message that Jordan Block cannot be built:

(a) n=0;

(b) n=-2;

(c) n=3.5;

(d) n=-2.5;

(e) n=4;

EXERCISE2 (5 points) Difficulty: Moderate

DESCRIPTION: In this part of the project, you will create a function added in a file. You will

code the sum of two matrices A and B, when they are of the same size, by using any of these

two definitions: (1) the sum of A and B is the matrix C whose columns are the sums of the

corresponding columns of A and B or (2) the sum of A and B is the matrix C whose entries

are the sums of the corresponding entries of the matrices A and B. Then, you will verify if

your output matrix C matches the output of a MATLAB built-in function for the sum, A+B.

**Create the function which begins with:

function C = added(A,B)

**First, your function has to verify whether the input matrices A and B have the same size. If

not, output a message ‘the matrices are not of the same size and cannot be

added’, and assign an empty matrix to C. After that, the program terminates.

**If the matrices can be added, you will continue with calculating the sum C of A and B using

one of the two definitions of the sum of two matrices (see above) – to code it, you can employ

“for” loop or vectorized statement. Output and display the matrix C.

**Next, you will use a logical “if” statement to verify whether the calculated matrix C

matches the output for a built-in MATLAB function A+B. If the outputs C and A+B do not

match, program an output message ‘check the code!’ – use disp or fprintf commands

– and, if you receive this message, you will need to correct the code and re-run Section.

This will be the end of your function added. Save it in your Current folder in MATLAB.

The rest of the work has to be done in the Live Script.

**Print your function in the Live Script:

type added

**Run the function added on the following sets of variables:

(a) A=magic(3), B=ones(4)

(b) A=ones(3,4), B=ones(3,3)

(c) A=randi(100,3,4), B=randi(100,3,4)

**You should receive a non-empty output matrix C in part (c), and you will continue working

in the Live Script with the matrices A, B, and C from part (c). You will “demonstrate” by

means of your function added and logical statements that the Commutative and Distributive

properties of the matrix addition hold (see (1)-(2) below):

3

(1) You will verify that the sum of A and B is the same as the sum of B and A. If it is the case,

the program should output a message that ‘commutative property holds for the

given A and B’.

(2) You will verify that the product of a scalar k by the matrix C, kC, is the same as the sum

of kA and kB. If this is the case, your program will output a message that ‘distributive

property holds for the given A and B’.

For this part you will need to input a scalar:

k=fix(10*rand(1))+5

Note: Make sure that you will receive the desired output messages for (1) and (2) above.

EXERCISE3 (5 points) Difficulty: Moderate

Theory: The transformation defined by the matrix

cos sin

sin cos

c s

A

s c

θ θ

θ θ

− − = =

performs a

rotation in 1 2 x x − plane through the angle θ radians in the counterclockwise direction for a

positive angleθ . The same type of a “plane” rotation in i j x x – plane can be performed in

( 2, 1 ) n n i jn ≥ ≤< ≤ , and it is called a Givens rotation.

A Givens rotation is represented by a matrix G which is defined by the following parameters:

n n( 2) ≥ , theta = θ (the angle of rotation), and i j , (1≤< ≤ i jn ). A matrix G can be created

by starting with the n n × identity matrix, that is G=eye(n), and, then, assigning the entries of

the matrix A to the entries of G at the intersections of rows i j , with columns i j , (i j < ) as

follows:

G(i,i)=c; G(i,j)=-s; G(j,i)=s; G(j,j)=c;

Geometrical Meaning: The matrix-vector product of G and a column vector x in n

performs a rotation of x in the i j x x -plane through θ radians, affecting, possibly, only rows

(=entries) i and j in x.

**Write a MATLAB function which begins with the line

function G=givensrot(n,i,j,theta)

and produces an n n × Givens matrix G under the conditions that 1≤< ≤ i jn and n ≥ 2 –

these conditions have to be verified in your code.

**If at least one of the conditions does not hold, the program outputs an empty matrix

G=[];

and terminates with an output message that Givens rotation matrix cannot be constructed.

**If all conditions hold, output and display the matrix G, as it specified in the Theory above.

**Print the function givensrot in your Live Script. type givensrot

**Run the function G=givensrot(n,i,j,theta) on each set of the variables in (1) – (5).

(Type pi for π in MATLAB.)

4

(1) n=1;i=1;j=2;theta=pi

(2) n=4;i=3;j=2;theta=pi/2

(3) n=5;i=2;j=4;theta=pi/4

(4) n=2;i=1;j=2;theta=-pi/2

(5) n=3;i=1;j=2;theta=pi

The rest of the work on this exercise has to be done in the Live Script using the output

matrix G from part (5):

**Based on the Geometrical Meaning of the product of G and a vector, predict what would be

the matrix (denote it GI) whose columns are the images of the vectors 123 eee , , , respectively,

under the Givens rotation defined in part (5); here, 123 eee , , are the columns of the 3 x 3

identity matrix I=eye(3), that is,

1

1

0

0

=

e , 2

0

1

0

=

e , 3

0

0

1

=

e .

Type and display GI in the Live Script.

**Then, calculate the actual products of G and each of the vectors 123 eee , , (in the order

indicated) – these products will form the columns of the matrix G*I (output and display G*I).

**Run a logical statement that would verify that your predicted matrix GI and the calculated

matrix G*I match within the margin 10^(-7). If the matrices match, output a message that your

prediction is correct – please make sure that you will receive this message.

Important: to compare two matrices, you need to use the function closetozeroroundoff

with p=7 (this function was created in Project 0).

**Print closetozeroroundoff in your Live Script.

type closetozeroroundoff

**Next, complete one more task. Input a vector:

x=ones(3,1);

Output and display the image of this vector under the transformation defined by the matrix G.

EXERCISE4 (5 points) Difficulty: Easy

We can define a Toeplitz matrix as an m n × matrix in which each descending diagonal from

left to right is constant. For example,

4321

5432

6543

A

=

is a 3 by 4 Toeplitz matrix.

Part 1. In this part, we will produce an m n × Toeplitz matrix using a vector a with

( 1) m n + − entries by assigning

Aij n i j i m j n ( , ) = +− = = a( ) 1: , 1: .

In the example above, a 3 4 × matrix A has been generated by the vector a = [1, 2,3, 4,5,6].

5

In this part of the exercise, you will write a code that creates an m n × Toeplitz matrix, when

possible. Then, you will run the function on the sets of variables m, n, and a to get the

corresponding Toeplitz matrices.

**Create a function in the file that begins with

function A=toeplitze(m,n,a)

Here, the inputs m and n will define the size of the Toeplitz matrix.

**First, your function has to check whether the number of entries in the vector a is

( 1) m n + − . If it is not the case, output a message that the dimensions mismatch and assign an

empty matrix to A. If vector a has exactly ( 1) m n + − entries, output and display the Toeplitz

matrix according to the formula given above.

This is the end of the function toeplitze.

**Print the function toeplitze in your Live Script. type toeplitze

**Run the function A=toeplitze(m,n,a)

on the following sets of variables (display the vector a for each of the below):

(a) m=4; n=2; a=1:5

(b) m=4; n=3; a=1:5

(c) m=4; n=3; a=1:7

(d) m=3; n=4; a=randi(10,1,6)

(e) m=4; n=4; a=[zeros(1,3), 1:4]

**Next, proceed in your Live Script with the following tasks:

(1) Type the vector a which defines a 6 by 6 upper triangular Toeplitz matrix with random

integer entries in the range between 0 and 100. Display the vector a. Run the function

A=toeplitze(m,n,a)on your a and the variables m,n indicated in this part.

(2) Output a 5 by 5 identity matrix by running the function A=toeplitze(m,n,a)on the

corresponding set of variables. Display the vector a (which you need to construct first).

Part 2. MATLAB has a built-in function, called toeplitz, which takes a different approach

to constructing a Toeplitz matrix. The specifications are on this page:

https://www.mathworks.com/help/matlab/ref/toeplitz.html?s_tid=doc_ta

**(a) Output a (symmetric) Toeplitz matrix T by running a MATLAB built-in function

T=toeplitz(r)

with

r=1:5;

**Run a logical statement in the Live Editor to verify that T is, indeed, a symmetric matrix.

Output the corresponding message if it is the case.

**(b) By running the function T=toeplitz(c,r), output a 6 by 5 Toeplitz matrix with the

same r and your choice of the vector c which should allow you to avoid the MATLAB

“warning” about “conflicting” entries.

6

More on Matrices; Application

EXERCISE5 (6 points) Difficulty: Hard

Theory: A vector with nonnegative entries is called a probability vector if the sum of its

entries is 1. A square matrix is called right stochastic if its rows are probability vectors; a

square matrix is called left stochastic if its columns are probability vectors; and a square

matrix is called doubly stochastic if both, the rows and the columns, are probability vectors.

**Write a MATLAB function that begins with

function [S1,S2,L,R]=stochastic(A)

L=[];

R=[];

It accepts as the input a square matrix A with nonnegative entries. Outputs L and R will be the

left stochastic and the right stochastic matrices generated, when possible, according to the

instructions given below. You will also calculate and display the vectors S1=sum(A,1)and

S2=sum(A,2) with the corresponding messages:

fprintf(‘the vector of sums down each column is\n’)

S1=sum(A,1)

fprintf(‘the vector of sums across each row is\n’)

S2=sum(A,2)

Then, you will use a conditional statement to proceed with the tasks outlined below. You may

also find it helpful to employ a logical command all.

**First, your function has to check whether A contains both a zero column and a zero row. If

yes, output a message “A is neither left nor right stochastic and cannot be scaled to either of

them”. The outputs L and R will be empty matrices (as assigned previously) – the empty

outputs should be suppressed.

**Then, the function checks whether A is: (1) doubly stochastic (assign: L=A;R=A); or (2)

only left stochastic (assign: L=A; R will stay empty), or (3) only right stochastic (assign R=A;

L will stay empty). In each of these cases, also output a message that comments on the type of

the matrix A.

**Finally, we consider a possibility that A is neither left nor right stochastic, but can be scaled

to the left stochastic and/or to the right stochastic. You will output a message “A is neither left

nor right stochastic but can be scaled to a stochastic matrix” and proceed with the scaling in

the ways outlined below:

(1) If neither S1 nor S2 has a zero entry, we are scaling* A to the left stochastic matrix L and

the right stochastic matrix R. You will also check if it is the case that the matrices L and R are

equal – in this case, A has been scaled to a doubly stochastic matrix.

(2) If S1 does not have a zero entry but S2 does, we scale A only to the left stochastic matrix

L (R stays empty).

(3) And, if S2 does not have a zero entry but S1 does, we scale A only to the right stochastic

matrix R (L stays empty).

Notes: In each of the cases, the non-empty outputs have to be displayed with the

corresponding messages. For the case of a doubly stochastic matrix, output the corresponding

message and display either L or R – to determine if L=R, you will need to use the function

closetozeroroundoff with p=7. ( This function was created in Project 0.)

7

*Scaling: To scale A to the left stochastic matrix L, we use vector S1 and multiply each

column of A by the reciprocal of the corresponding entry of S1. To scale A to the right

stochastic matrix R, we use the vector S2 and multiply each row of A by the reciprocal of the

corresponding entry of S2.

**Print the functions stochastic, jord, and closetozeroroundoff in your Live Script.

(The function jord was created in Exercise 1 of this project.)

**Run the function [S1,S2,L,R]=stochastic(A) on each of the matrices below (display

the input matrices in your Live Script):

(a) A=[0.5,0,0.5,0; 0,0,1,0;0.5,0,0.5,0;0,0,0,1]

(b)A = transpose(A)

(c)A=[0.5, 0, 0.5; 0, 0, 1; 0, 0, 0.5]

(d)A=transpose(A)

(e)A=[0.5, 0, 0.5; 0, 0.5, 0.5; 0.5, 0.5, 0]

(f)A=magic(3)

(g)B=[1 2;3 4;5 6]; A=B*B’

(h)A=jord(5,4)

(k)A=randi(10,5,5);A(:,1)=0;A(1,:)=0

NOTE: Please make sure you will verify that all your outputs and messages match the

corresponding definitions of the stochastic matrices. If not, make corrections in your code!

EXERCISE6 (5 points) Difficulty: Moderate

In this exercise, you will use Newton’s method to approximate a real zero of a given

function.

Theory: A number x is a zero of a function f if f x( ) = 0 . A real zero is an x-intercept of the

function. Also, the zeros of a function f are the roots of the equation f x( ) = 0 .

To approximate a real zero of the function f x( ) by Newton’s method, we, first, choose an

initial approximation 0 x of the specified zero. The consecutive N iterations 1 2 , ,…, N xx x are

defined by the formula:

( )

( ) 1 , 0: . ‘

n

n n

n

f x xx nN

f x + =− =

Note: if f x ‘( ) is close to zero at an initial point 0 x , the process will converge very slowly. It

makes Newton’s method sensitive to the choice of an initial approximation.

In this exercise, you will work with two functions:

Fx x x ( ) = +− arctan( ) 1 and ( ) 3 Gx x x = −−1.

**First, you will create a function handle. Please refer to the documentation in MATLAB:

https://www.mathworks.com/help/matlab/matlab_prog/creating-a-functionhandle.html?searchHighlight=handle&s_tid=doc_srchtitle#buvu9u8-1

We will create handles to anonymous functions and calculate derivatives of the function

handles as it follows.

8

**Begin with typing in the Live Script:

format

format compact

syms x

F = @(x) atan(x) + x – 1

F1 = eval([‘@(x)’ char(diff(F(x)))])

G=@(x) x.^3-x-1

G1=eval([‘@(x)’ char(diff(G(x)))])

Notes: (1) a MATLAB command syms x defines a symbolic variable x;

(2) F and G are the function handles and F1 and G1 are the function handles for the first

derivatives of F and G, respectively.

**Next, we create and output 2-D plots of the functions F, G, and y=0 on the interval [-2, 2]

in order to visualize the x-intercepts and choose initial approximations. An initial value 0 x

should be chosen close to the x-intercept which we are approximating.

**Type in the Live Script the code given below – it will output the graphs of F and G together

with the function y=0 after you Run the section.

yzero=@(x) 0.*x.^(0)

x=linspace(-2,2);

plot(x,F(x),x,yzero(x));

plot(x,G(x),x,yzero(x));

Note: here we have created another symbolic function yzero.

It is obvious from the properties of the function F(x) that it has only one real zero, which we

will approximate. Concerning the function G, which is a polynomial of the third degree (we

will denote it p), we are going to verify that it has only one real zero. In order to do that, we

find all zeros of the polynomial p using MATLAB built-in functions sym2poly and roots:

the function sym2poly(p) outputs the vector of the coefficients of the polynomial p (in

descending order according to the degree), and the composition of two functions

roots(sym2poly(p))outputs all zeros of the polynomial p.

**Type in the Live Script:

syms x

p=x^3-x-1;

roots(sym2poly(p))

After you Run Section, this part of the code will output the three zeros of the polynomial p –

two of them are complex conjugate numbers and one is a real zero that we will approximate.

Next, we proceed with constructing a function in the file that approximates a real zero.

**Create a function called newtons. It begins with:

function root=newtons(fun,dfun,x0)

format long

9

The inputs fun and dfun are the function and its first derivative, respectively, and x0 is the

initial approximation. The output root will be our approximation of the real zero of a

function. We will program consecutive iterations according to the Newton’s method, and we

will assign to root the iteration which will be the first one falling within a margin of 10^(-12)

from the MATLAB approximation x of that zero – the last is delivered by a built-in

MATLAB function fzero. The details are below:

**Type the line

x=fzero(fun,x0)

in your function newtons and output x with a message that it is a MATLAB approximation of

the real zero of the function.

**Then, your function newtons will calculate consecutive iterations , 0: n xn N = , using

Newton’s Method (see Theory above). You can employ a “while loop” here. The loop will

terminate when, for the first time, abs( Nx -x)<10^(-12) for some consecutive iteration Nx .

Output with a corresponding message the number of iterations N, and assign the last iteration,

Nx , to the output root.

This will be the end of your function newtons.

**Print the function newtons in your Live Script.

**Next, proceed with the following tasks in the Live Script:

Part (a)

You will work with the function F in this part. Using the graph of the function F, choose three

different values of the initial approximation x0 of the zero of F.

**Then, input the function handles:

fun=F;

dfun=F1;

and run the function root=newtons(fun,dfun,x0) for each of your three choices of x0 –

one at a time.

Part (b)

You will work with the function G in this part.

**First, input the corresponding function handles.

**Then, run root=newtons(fun,dfun,x0) for each of the initial approximations (1)-(8) of

the real zero of G:

(1) x0=1.3;

(2) x0=1;

(3) x0=0.6;

(4) x0=0.577351;

(5) Pick the initial value x0 to be the positive zero of the derivative of G(x). Type

x0=1/sqrt(3)

(display x0) and run the function root=newtons(fun,dfun,x0).

(6) x0=0.577;

(7) x0=0.4;

(8) x0=0.1;

See next page for the BONUS on this Exercise.

10

Analyze the patterns of the consecutive iterations after you run your function newtons on the

choices of the initial approximations (1)-(8) in Part (b).

BONUS: 1 point!

% Describe the general pattern for the numbers of iterations which you observed for the

choices (1)-(8), excluding choices (4)-(6) out of this sequence.

BONUS: 1 point!

% What did you notice when you run your function on choices (4)-(6)? How would you

explain that “strange” pattern for the consecutive iterations?

Hint: to get the bonus points, you could temporarily display all intermediate iteration that

come from the “loop”. You may also need to review the Newton’s Method theory in more

detail.

THIS IS THE END OF PROJECT 1