I am trying to write a function in Julia that solves the nonlinear equation f(x)=0 using Newton iteration. I am very much a beginner in Julia so bear with me. In the assignment, my instructor provided this first line of code:
newton(f, x_0; maxiterations=10, tolerance=1e-14, epsilon=1e-7)
He also provided this statement: "the optional epsilon argument provides the value of ε used in the finite-difference approximation f′(x) ≈ (f(x + ε) − f(x)) / ε." In the past, I've created a function for Newton's Method in MATLAB but I assumed that the first derivative of f(x) must be one of the function's inputs. In this assignment, it appears that he wants me to use this approximation formula. Anyways, here is my code so far.
function newton(f,x_0; maxiterations=10, tolerance=1e-14, epsilon=1e-7)
x = [x_0] # assign x_0 to x and make x a vector
fd = (f(x.+epsilon) - f(x))./epsilon # fd is the first derivative of f(x), calculated from
# the finite-difference approximation
# create for loop to begin iteration
for n = 0:maxiterations
if abs(x[n+1]-x[n]) < tolerance # if the absolute value of the difference of x[n+1] and x[n]
# is less than the tolerance, then the value of x is returned
return x
end
if abs(f(x[n])) < tolerance # if the absolute value of f(x[n]) is less than the tolerance,
# then the value of x is returned
return x
end
push!(x, x[n] - (f(x[n]))/fd) # push each calculated value to the end of vector x,
# and continue iterating
end
return x # after iteration is complete, return the vector x
end
After executing this function, I defined the equation which should be used to determine sqrt(13) and called the newton function with an initial guess of x_0=3.
f(x) = x^2 - 13
newton(f,3)
Here is the error message I'm encountering after the newton function is called:
MethodError: no method matching ^(::Vector{Float64}, ::Int64)
Closest candidates are:
^(::Union{AbstractChar, AbstractString}, ::Integer) at strings/basic.jl:730
^(::LinearAlgebra.Hermitian, ::Integer) at /Applications/Julia-1.8.app/Contents/Resources/julia/share/julia/stdlib/v1.8/LinearAlgebra/src/symmetric.jl:696
^(::LinearAlgebra.Hermitian{T, S} where S<:(AbstractMatrix{<:T}), ::Real) where T at /Applications/Julia-1.8.app/Contents/Resources/julia/share/julia/stdlib/v1.8/LinearAlgebra/src/symmetric.jl:707
...
Stacktrace:
[1] literal_pow
@ ./intfuncs.jl:340 [inlined]
[2] f(x::Vector{Float64})
@ Main ./In[35]:3
[3] newton(f::typeof(f), x_0::Int64; maxiterations::Int64, tolerance::Float64, epsilon::Float64)
@ Main ./In[34]:5
[4] newton(f::Function, x_0::Int64)
@ Main ./In[34]:2
[5] top-level scope
@ In[35]:4
[6] eval
@ ./boot.jl:368 [inlined]
[7] include_string(mapexpr::typeof(REPL.softscope), mod::Module, code::String, filename::String)
@ Base ./loading.jl:1428
In this problem, I'm supposed to make my function return a vector of successive approximations and make a logarithmic error plot that includes a curve for the theoretical error estimate for Newton iteration. I would appreciate any guidance in correcting the issue(s) in my code so I can move on to create the error plot.
Thank you.
You can use
newton(x->x.^2 .- 13, 3.)
to manipulate vectors instead of scalars. But your newton code is still bugged!
you start with zeros indexes(not good) and x[n+1] is used before to be allocated and it lacks some mathematical steps. I can propose you
function newton(f,x_0; maxiterations=10, tolerance=1e-14, epsilon=1e-7)
x = x_0
iterates = [x]; # iterates will store the iterates' sequence
# create for loop to begin iteration
for n =1:maxiterations
fd = (f(x.+epsilon) - f(x))./epsilon
tmp=x - f(x) / fd
if (abs(tmp-x) < tolerance || abs(f(x)) < tolerance )
break
end
x = tmp
push!(iterates,x)
end
return iterates
end
Now the result seems to be good
julia> newton(x->x.^2 .- 13, 3.)
5-element Vector{Float64}:
3.0
3.6666666569022053
3.6060606066709835
3.605551311439882
3.60555127546399