Hello, dear friend, you can consult us at any time if you have any questions, add WeChat: daixieit

MATH50003 Numerical Analysis (2021-2022) Practice Computer-based Exam

For each problem, replace the # TODO to complete the question. The unit tests are provided to help you test your answers. You have 1 hour to complete the exam, as well as 1 hour for downloading/uploading.

Problems are marked A/B/C to indicate difficulty ("A" being most difficult). Partial credit will be awarded for reasonable attempts even if the tests are not passed.

You may use existing code from anywhere but you are REQUIRED to cite the source if it is not part of the module material, ideally by including a weblink in a comment You MUST NOT ask for help online or communicate with others within or outside the module. Failure to follow these rules will be considered misconduct

You should use the following packages:

In [1]:

using LinearAlgebra, SetRounding, Test

WARNING It may be necessary to restart the kernel if issues arise. Remember to reload the packages when you do so.

1. Numbers

Problem 1.1 (C) Complete the following function divideby3 (x) that returns a tuple a, b such that a is the largest Float64 less than or equal to x/3 and b is the smallest Float64 greater than or equal to x/3 .

In [2]:

function divideby3(x)

# TODO: assign a, b so that a W x W b where b is either equal to a or the next float

a, b

end

Out[2]:

divideby3 (generic function with 1 method)


In ⑶:

x 0.1 # arbitary x

a, b = divideby3(x)

©test a W big(x)/3 W b

©test b == a I I b == next float (a)

UndefVarError: a not defined

Stacktrace:

[1] divideby3(x::Float64)

@ Main . /In[2]:4

[2] top-level scope

@ In[3]:2

[3] eval

@ . /boot.jl:373 [inlined]

[4] include_string(mapexpr::typeof(REPL. softscope), mod::Module, code::String, file name::String)

@ Base . /loading, jl:1196

2. Differentiation

Problem 2.1 (C) Use the following off-center finite-difference approximation a gm

J 3h

with an appropriately chosen h to approximate f(x) = cos(x2)

at x = 0.1 to 5 digits accuracy.

In [4]:

function fd(x)

# TODO: implement a finite-difference rule

# to approximate f' (x)

# for f (x) = cos(x"2)

# with step-size h chosen to get sufficient accuracy

end

Out [4]:

fd (generic function with 1 method)


In [5]:

©test abs(fd(0. 1) + 2*0. l*sin(0. 12)) W IE-5

Error During Test at In[5] : 1

Test threw exception

Expression: abs (fd(O. 1) + 2 * 0. 1 * sin(0. 1 “ 2)) W 1. 0e~5

MethodError: no method matching +(::Nothing, ::Float64)

Closest candidates are:

+ (::Any, ::Any, ::Any, ::Any...) at /Applications/Julia-1. 7. app/Contents/Resourc es/julia/share/julia/base/operators, jl:655

+ (::Base. TwicePrecision, ::Number) at /Applications/Julia-1. 7. app/Contents/Resou rces/julia/share/julia/base/twiceprecision. jl:279

+ (::UniformScaling, ::Number) at /Applications/Julia-1. 7. app/Contents/Resources/ julia/share/julia/stdlib/vl. 7/LinearA1gebra/src/uniformsealing. j1:145

Stacktrace:

[1] top-level scope

@ /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/vl. 7/

Test/src/Test. _il:445

[2] eval

@ ・/boot, i1:373 [inlined]

[3] include_string (mapexpr::typeof(REPL. softscope), mod::Module, code::String, filename::String)

@ Base . /loading, jl:1196

[4] sof tscope_include_string (m::Module, code::String, filename::String)

© SoftGlobalScope julia/packages/SoftGlobalScope/uWzH/src/SoftGlobalScope^j. 1:65

[5] execute_request (socket::ZMQ. Socket, msg:: IJulia. Msg)

@ IJulia "7. iulia/Dackases/Llulia/e8kciU/src/execute request. jl:67

[6] #invokelatest#2

@ . /essentials, jl:716 [inlined]

[7] invokelatest

@ . /essentials, jl:714 [inlined]

[8] event loop (socket:: ZMQ. Socket)

@ IJulia ^/. jul i a/packages/I.Tul ia/ e8kqU/src/event loop, .i 1:8

[9] (:: IJulia. var"#15#18z/)()

@ IJulia /task, jl:423

There was an error during testing

Stacktrace:

[1] record(ts::Test. FallbackTestSet, t::Union{Test. Error, Test. Fail})

@ Test /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/v

1. 7/Test/src/Test. jl:903

[2] do_test(result::Test. ExecutionResult, orig_expr::Any)

@ Test /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/v 1.7/Test/src/Test. jl:637

[3] top-level scope

@ /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/vl. 7/Te st/src/Test. jl:445

[4] eval

@ . /boot, jl:373 [inlined]

[5] include_string(mapexpr::typeof(REPL. softscope), mod::Module, code::String, file name::String)

@ Base . /loading, jl:1196

Problem 2.2 (A) Consider a 2D version of a dual number:

a + bex + cey

such that

€?=€§= Cx 与=0・

Complete the following implementation supporting + and * (and assuming a, b, c are Float64 ). Hint: you may need to work out on paper how to multiply (s. a + s. b e_x + s. c €_y)*(t. a + t. b e_x + t. c €_y) using the relationship above.

In [6]:

import Base: *, +, struct Dual2D

a

b

c

end


function + (s: :Dual2D, ## TODO: Implement end

t::Dual2D)

+, returning a Dual2D


s::Dual2D)

c * Dual2D(.returning a Dual2D

end

t::Dual2D)

Dual2D(...) * Dual2D(...), returning a Dual2D

en

Out[6]

* (generic function with 366 methods)


In [7]:

f function (x, y) # (x+2y"2)"3 using only * and +

z = (x + 2y * y)

z * z * z

end

x, y = 1., 2.

©test f (Dual2D(x, 1., 0.), Dual2D(y, 0., 1. )) == Dual2D (f (x, y), 3(x+2y"2)”2, 12y* (x+2y"2) ”2)

# This has computed the gradient as f (x, y) + f_x*€Lx + f_y*€Ly

# == (x+2y"2)"3 + 3(x+2y"2)"2*e_x + 12y(x+2y"2)"2*6_y

Error During Test at In [7] : 7

Test threw exception Expression: f (Dual2D(x, 1. 0, 0. 0), Dual2D(y, 0. 0, 1. 0)) = Dual2D(f(x, y), 3 * (x

+ 2 * y " 2) "2, (12y) * (x + 2 * 广 2) "2)

MethodError: no method matching *(::Nothing, ::Dual2D)

Closest candidates are:

* (::Any, ::Any, ::Any, ::Any...) at /Applications/Julia-1. 7. app/Contents/Resourc es/julia/share/julia/base/operators, j1:655

*(::Number, ::Dual2D) at In[6]:13

*(::Dual2D, ::Dual2D) at In[6]:18

Stacktrace:

[1] (: : var"#l#2”) (x:: Dual2D, y: : Dual 2D)

@ Main ./In[7]:2

[2] top-level scope

@ /App1ications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/vl. 7/Test/src/Test. .11:445

[3] eval

@ . /boot. .11:373 [inlined]

[4] include_str ing (mapexpr::typeof(REPL. softscope), mod::Module, code::String, filename::String)

@ Base ./loading, jl:1196

[5] sof tscope_include_string (m::Module, code::String, filename::String)

@ SoftGlobalScope ~/. julia/Dackages/SoftGlobalScoDe/u4UzH/src/SoftGlobalScope.

J1165

[6] execute_request (socket::ZMQ. Socket, msg:: IJulia. Msg)

@ IJulia "7. iulia/Dack&zes/I.Tulia/e8kaU/src/execute request, jl:67

[7] #invokelatest#2

@ ./essentials. il:716 [inlined]

[8] invokelatest

@ ・/essentials, jl:714 [inlined]

[9] event loop (socket:: ZMQ. Socket)

@ IJulia V. julia/packages/IJulia/e8kqU/src/eventlooph_jjj_8

[10] (:: IJulia. var"#15#18")()

@ IJulia . /task, il:423

There was an error during testing

Stacktrace:

[1] record(ts::Test. FalIbackTestSet, t::Union{Test. Error, Test. Fail})

@ Test /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/

vl. 7/Test/src/Test. jl:903

[2] do_test(result::Test. ExecutionResult, orig_expr::Any)

@ Test /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/

vl. 7/Test/src/Test. jl:637

[3] top-level scope

@ /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/vl. 7/

Test/src/Test. j1:445

[4] eval

© . /boot, jl:373 [inlined]

[5] include_string(mapexpr::typeof(REPL. softscope), mod::Module, code::String, fi 1ename::String)

@ Base . /loading, jl:1196

3. Structured Matrices

Problem 3.1 (C) Add an implementation of inv (:: PermutationMatrix) to return the inverse permutation as a PermutationMatrix . Hint: use invperm .

In [8]:

import Base: get index, size, *, inv

struct PermutationMatrix < : AbstractMatrix (Int}

p::Vector{Int) # represents the permutation whose action is v[p] function PermutationMatrix(p::Vector)

sort (p) = 1: length(p) | | error (zzinput is not a valid permutation") new (p)

end

end

size(P::PermutationMatrix) = (length(P. p), length(P. p))

getindex(P::PermutationMatrix, k::Int, j::Int) = P. p[k] = j ? 1 : 0

*(P::PermutationMatrix, x::AbstractVector) = x[P. p]

function inv(P::PermutationMatrix)

# TODO: return a PermutationMatrix representing the inverse permutation

end

Out ⑻:

inv (generic function with 33 methods)


In [9]:

P PermutationMatrix([3, 4, 2, 1])

©test inv(P) isa PermutationMatrix

©test P* inv (P) == I

Test Failed at In[9]:2

Expression: inv(P) isa PermutationMatrix

Evaluated: nothing isa PermutationMatrix

There was an error during testing

Stacktrace:

[1] record(ts::Test. FalIbackTestSet, t::Union{Test. Error, Test. Fail))

@ Test /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/v 1.7/Test/src/Test. jl:903

[2] do_test(result::Test. ExecutionResult, orig_expr::Any)

@ Test /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/v 1.7/Test/src/Test. jl:637

[3] top-level scope

@ /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/vl. 7/Te st/src/Test. jl:445

[4] eval

@ . /boot, jl:373 [inlined]

[5] include_string(mapexpr::typeof(REPL. softscope), mod::Module, code::String, file name::String)

@ Base . /loading, jl:1196

4. Decompositions

Problem 4,1 (C) Forx G 2 implement the reflection defined by

y = x+ ||x||e„
w = y/||y||

Qx \= I 2wwT

in lowerhouseholderreflection(x) , which should return a Matrix (Float64} . You may assume that x is a Vector(Float64).

In [10]:

function 1owerhous eho1derr e f1e c t i on(x)

## TODO: implement the householder reflector defined above

end

Out[10]:

lowerhouseholderreflection (generic function with 1 method)

In [11]:

x = [1.0, 2, 3, 4]

Q = lowerhouseholderreflection(x) ©test Q*x 2 [zeros(3); -norm(x)] ©test Q' Q Q I

©test Q e Q

Error During Test at In [11] : 3

Test threw exception

Expression: Q * x R [zeros(3); -(norm(x))]

MethodError: no method matching *(::Nothing, ::Vector{Float64))

Closest candidates are:

*(::Any, ::Any, ::Any, ::Any...) at /Applications/Julia-1. 7. app/Contents/Resourc es/julia/share/julia/base/operators. j1:655

*(::StridedMatrix{T}, ::StridedVector{S}) where {T<:Union(Float32, Float64, Comp lexF32, ComplexF64), S:Real} at /Applications/Julia-1. 7. app/Contents/Resources/juli a/share/julia/stdlib/vl. 7/LinearA1gebra/src/matmu1. jl:44

*(::Union(SparseArrays. AbstractSparseMatrixCSC{TA, Ti}, SubArray{TA, 2, <:Sparse Arrays. AbstractSparseMatrixCSC{TA, Ti}, Tuple{Base. Slice(Base. 0neTo{Int64}), 1}} whe re I<:AbstractUnitRange} where Ti, ::Union(StridedVector, BitVector}) where TA at /A ppiications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/vl. 7/SparseArr ays/src/linalg. jl:50

Stacktrace:

[1] top-level scope

@ /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/vl. 7/ Test/src/Test. il:445

[2] eval

@ . /boot. il:373 [inlined]

[3] include_string (mapexpr::typeof(REPL. softscope), mod::Module, code::String, filename::String)

@ Base ./loading, jl:1196

[4] sof tscope_include_string (m::Module, code::String, filename::String)

@ SoftGlobalScope ~/. iulia/Dackages/SoftGlobalScoDe/u4UzH/src/SoftGlobalScope. j 1:65

[5] execute_request (socket::ZMQ. Socket, msg:: IJulia. Msg)

@ IJulia ~/. julia/Dackages/LTulia/e8kaU/src/execute request, jl:67

[6] #invokelatest#2

@ ./essentials. il:716 [inlined]

[7] invokelatest

@ ./essentials, jl:714 [inlined]

[8] event loop (socket:: ZMQ. Socket)

@ IJulia ~/, julia/packages/I.Tulia/e8kqU/src/eventloop, jl:8

[9] (:: IJulia. var"#15#18")()

@ IJulia . /task. 11:423

There was an error during testing

Stacktrace:

[1] record(ts::Test. FalIbackTestSet, t::Union{Test. Error, Test. Fail})

@ Test /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/v 1. 7/Test/src/Test. jl:903

[2] do_test(result::Test. ExecutionResult, orig_expr::Any)

@ Test /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/v 1.7/Test/src/Test. jl:637

[3] top-level scope

@ /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/vl. 7/Te st/src/Test. jl:445

[4] eval

@ . /boot, jl:373 [inlined]

[5] include_string(mapexpr::typeof(REPL. softscope), mod::Module, code::String, file name::String)

@ Base . /loading, jl:1196

Problem 4.2 (A) Complete the function ql (A) that returns a QL decomposition, that is5 Q is an orthogonal matrix and L is lower triangular. You may assume that A is a square Matrix {Float64) . Hint: use Problem 4.1 to lower triangularise.

In [12]:

function ql(A)

m, n size (A) m == n I I error ("not square")

## TODO Create Q and L such that Q'Q = I and L is lower triangular

end

Out [12]:

ql (generic function with 1 method)

In [13]:

A =[1.0 2 3; 1 4 9; 1 1 1]

Q, L ql(A)

©test Q' Q e I

@test Q*L 2 A

©test L 2 tril(L) # it is acceptable to have small non-zero entries in L

BoundsError: attempt to access Bool at index [2]

Stacktrace:

[1] indexediterate(I::Bool, i::Int64, state::Nothing)

@ Base ・/tuple. jl:98

[2] top-level scope

@ In[13]:2

[3] eval

© . /boot, jl:373 [inlined]

[4] include_string(mapexpr::typeof(REPL. softscope), mod::Module, code::String, file name::String)

@ Base . /loading, jl:1196

5. Singular Value Decomposition

Problem 5.1 (C) Find the best rank-4 approximation (in the 2-norm) to

/(x, y) = cos(x y)/(x + y+l) evaluated at an evenly spaced 100 x 100 grid on the square [0,1]2.

In [14]:

function bestrank4()

# TODO: return best rank-4 approximation

end

Fr = bestrank4 ()


In [15]:

x  9/99

y = 10/99

©test rank(Fr) == 4

©test abs (Fr [10, 11] - cos (x - y)/ (x + y + 1)) W 2E-6

Error During Test at In[15] : 3

Test threw exception

Expression: rank(Fr) = 4

MethodError: no method matching rank(::Nothing)

Closest candidates are:

rank(::Number) at /Applications/Julia-1. 7. app/Contents/Resources/julia/share/jul ia/stdlib/vl. 7/LinearAlgebra/src/generic.jl:1013

rank (::CholeskyPivoted) at /Applications/Julia-1. 7. app/Cont ent s/Re sour ces/julia/

share/julia/stdlib/vl. 7/LinearAlgebra/src/cholesky. jl:697

rank(::SuiteSparse. SPQR. QRSparse) at /App1ications/Ju1ia-1. 7. app/Contents/Resour ces/julia/share/julia/stdlib/vl. 7/SuiteSparse/src/spqr. jl:324

Stacktrace:

[1] top-level scope

@ /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/vl. 7/

Test/src/Test. jl:445

[2] eval

@ ・/boot, jl:373 [inlined]

[3] include_string (mapexpr::typeof(REPL. softscope), mod::Module, code::String, filename::String)

@ Base . /loading, jl:1196

[4] sof tscope_include_string (m: :Module, code::String, filename::String)

@ SoftGlobalScope ~/. iulia/Dackages/SoftGlobalScoDe/u4UzH/src/SoftGlobalScope. i

1:65

[5] execute_request (socket::ZMQ. Socket, msg:: IJulia. Msg)

@ IJulia ~/・ julia/Dackages/L]ulia/e8kciU/src/execute request. jl:67

[6] #invokelatest#2

@ . /essentials, jl:716 [inlined]

[7] invokelatest

@ ・/essentials, jl:714 [inlined]

[8] event loop (socket:: ZMQ. Socket)

@ IJulia "/. iulia/packages/LTulia/eSkaU/src/eventlooD. _jl:8

[9] (:: IJulia. var"#15#18z)()

@ IJulia . /task, il:423

There was an error during testing

Stacktrace:

[1] record(ts::Test. FalIbackTestSet, t::Union{Test. Error, Test. Fail})

@ Test /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/v

1. 7/Test/src/Test. jl:903

[2] do_test(result::Test. ExecutionResult, orig_expr::Any)

@ Test /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/v

1. 7/Test/src/Test. jl:637

[3] top-level scope

@ /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/vl. 7/Te st/src/Test. jl:445

[4] eval

@ . /boot, jl:373 [inlined]

[5] include_string(mapexpr::typeof(REPL. softscope), mod::Module, code::String, file name::String)

@ Base . /loading, jl:1196

6. Differential Equations

Problem 6J (A) Complete the function airyai (n) that returns a length-n Vector (Float64}

U\

such that Uk approximates the solution to the equation

u(0) = 1

w(l) = 0
/ — xu = 0

at the point = (k — l)/(n 1) using the second order finite-difference approximation:

Uff(xk) « Uk_i-2"k+"k+l

for = 2,..., n - 1, in O(n) operations.

In [16]:

function airy(n)

# T0D0: return a Vector(Float64} approximating the solution to the ODE end

Out[16]:

airy (generic function with 1 method)


In [17]:

u airy(1000)

©test u[l] == 1

©test u[end] == 0

# this compares agianst the exact formula ©test abs(u[500] - 0.4757167332829094) W 2E-8

Error During Test at In[17] : 2

Test threw exception

Expression: u[l] = 1

MethodError: no method matching getindex(::Nothing, ::Int64)

Stacktrace:

[1] top-level scope

@ /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/vl. 7/ Test/src/Test. .11:445

[2] eval

@ ./boot, i1:373 [inlined]

[3] include_string (mapexpr::typeof(REPL. softscope), mod::Module, code::String, filename::String)

@ Base . /loading, il:1196

[4] sof tscope_include_string (m::Module, code::String, filename::String)

@ SoftGlobalScope ~/・ julia/packages/SoftGlobalScope/u4UzH/src/SoftGlobalScope, i 1:65

[5] exe cut e_re quest (socket:: ZMQ. Socket, msg:: I Julia. Msg)

@ IJulia /. julia/Dackases/I.lulia/eBkaU/src/execute request. jl:67

[6] #invokelatest#2

@ . /essentials, jl:716 [inlined]

[7] invokelatest

@ ・/essentials, il:714 [inlined]

[8] event loop (socket:: ZMQ. Socket)

@ IJulia ~/. julia/packages/IJulia/e8kqU/src/eventlooiL_jJd_8

[9] (: : IJulia. var"#15#18")() @ IJulia ./task, il:423

There was an error during testing

Stacktrace:

[1] record(ts::Test. FalIbackTestSet, t::Union(Test. Error, Test. Fail})

@ Test /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/v 1.7/Test/src/Test. jl:903

[2] do_test(result::Test. ExecutionResult, orig_expr::Any)

@ Test /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/v 1. 7/Test/src/Test. jl:637

[3] top-level scope

@ /Applications/Julia-1. 7. app/Contents/Resources/julia/share/julia/stdlib/vl. 7/Te st/src/Test. jl:445

[4] eval

@ . /boot.jl:373 [inlined]

[5] include_string(mapexpr::typeof(REPL. softscope), mod::Module, code::String, file name::String)

@ Base . /loading, jl:1196