MiniProjectCh8.sagews

# Mini project: Comparing convergence of the simple iteration method and the Newton's method
#Consider the function f(x)=exp(-x)*(x^2-5*x-2)+1.
# Plot the function to convince yourself that the function has two roots. For each of the roots use the "find_root" function to obtain high accuracy value of the root1 and root2.
#For each of the roots Do the following
#(a) (challenging!) Convert the root-finding problem for f to a fixed point problem g(x)=x with appropriate function g.
#(b) Choose an initial approximation of the root that guarantees convergence of the simple iteration and make a list L with n=7 iterates using simple iteration.
#(c) Use the same initial guess for solving the root-finding problem by Newton's method. Make a list N with n=7 iterates using Newton's iterations.
#(d) Make a table with 7 rows for x[n]-values, n=0..6, and 3 columns with the table header n, L[n],and N[n]. Compare convergence of these two methods.

reset

<function reset at 0x7effb05d00c8>

f(x)=exp(-x)*(x^2-5*x-2)+1
plot(f(x),x,-1,3, figsize=[3,4])

root1=find_root(f,-1,0);root1

-0.23054999321338798

root2=find_root(f,1,3);root2

2.0895848494844653

#Version 1
# (a) Solve the equation f(x)=0 for x. Hint: there are three options here. Do not solve for x in the exponential function. The other two options will work.
g1(x)=(x^2+exp(x)-2)/5

line_plot1=plot(x,-1,0)
g1plot=plot(g1,-1,3,figsize=[4,3])
g1prime_plot=plot(((2/5)*x+exp(x)/5),-1,3,color='red')
(g1plot+g1prime_plot+line_plot1).show()

#(b) g1 works well for root1; gprime near the root1 is ?. Let x0=-0.9
def simple_iteration(g,a,n):
    L=[a]
    temp=a
    for j in range(1,n+1):
        temp=g(L[j-1]).n()# Encode the iteration formula here
        L.append(temp)
    return L

simple_it_list1=simple_iteration(g1,-0.9,7);simple_it_list1

[-0.900000000000000, -0.156686068051880, -0.224095411245832, -0.230108474460778, -0.230520533881549, -0.230548030998561, -0.230549862530082, -0.230549984509961]

# (c) Plot of f shows that the same initial guess will work for approximating root1 by Newton's method ( since ?>0)

def my_newton(f,a,n):
    L=[a]
    fprime=diff(f,x)
    for j in range(n):
        xold=N(L[j])
        #Your code for Newton's formula is incorrect
        #xnew=f(L[j-1])-(f(L[j-1])/fprime(L[j-1])) # Encode Newton's formula here
        #Here is a correct code:
        xnew=xold-(f(x=xold)/fprime(x=xold))
        L.append(xnew)
    return L

newtons_list1=my_newton(f,-0.9,7);newtons_list1

[-0.900000000000000, -0.532386779451968, -0.313764802062772, -0.238603685771938, -0.230633323982670, -0.230550002235008, -0.230549993213389, -0.230549993213389]

#(d) Recall that the exact root1 is -0.23054999321338895
nlist=[j for j in range(7)]
table(columns=[nlist,simple_it_list1,newtons_list1],header_row=['n',"Simple iterations", "Newton's iterations"], frame=True)#learn the synatax of this valuable command

+---+--------------------+---------------------+ | n | Simple iterations | Newton's iterations | +===+====================+=====================+ | 0 | -0.900000000000000 | -0.900000000000000 | +---+--------------------+---------------------+ | 1 | -0.156686068051880 | -0.532386779451968 | +---+--------------------+---------------------+ | 2 | -0.224095411245832 | -0.313764802062772 | +---+--------------------+---------------------+ | 3 | -0.230108474460778 | -0.238603685771938 | +---+--------------------+---------------------+ | 4 | -0.230520533881549 | -0.230633323982670 | +---+--------------------+---------------------+ | 5 | -0.230548030998561 | -0.230550002235008 | +---+--------------------+---------------------+ | 6 | -0.230549862530082 | -0.230549993213389 | +---+--------------------+---------------------+

# Comment:the simple iterations its approching -.23and the newtons its approching -.23

###################################

g2(x)=(5*x+2-exp(x))^1/2

diff(g2(x))

-1/2*e^x + 5/2

# Version 2: Root 2
# (a)
line_plot2=plot(x,1,2.5)
g2(x)=(5*x+2-exp(x))^0.5
g2plot=plot(g2(x),1,2.5)
g2prime_plot=plot(-1/2*exp(x) + 5/2,1,2.5,color='red')# g2prime is incorrect
g2prime=diff(g2,x);g2prime
(g2plot+g2prime_plot+line_plot2).show(figsize=[4,3])

x |--> -0.500000000000000*(e^x - 5)/sqrt(5*x - e^x + 2)

# (b) g2 works well for root2; there is an interval containing root2 where gprime is ?. Let x0=1.2
simple_it_list2=simple_iteration(g2,1.2,7);simple_it_list2

[1.20000000000000, 2.16330374133256, 2.02896192145109, 2.13040446092206, 2.05760815511591, 2.11206412847833, 2.07244694614844, 2.10190300751602]

#c)
f(x)=(-1/2*exp(x) + 5/2)
newtons_list2=my_newton(f,1.2,7);newtons_list2

[1.20000000000000, 1.70597105956101, 1.61395085974027, 1.60944808047909, 1.60943791248579, 1.60943791243410, 1.60943791243410, 1.60943791243410]

#(d) Recal exact root2: 2.0895848494844653
table(columns=[nlist,simple_it_list2,newtons_list2],header_row=['n',"Simple iterations", "Newton's iterations"], frame=True)

###############
#Grading: There are 20 red question marks; correct replacement for each is 0.5 point. There are two brown question marks (theoretical); each replacement 2 points. Total: 12 points.

+---+-------------------+---------------------+ | n | Simple iterations | Newton's iterations | +===+===================+=====================+ | 0 | 1.20000000000000 | 1.20000000000000 | +---+-------------------+---------------------+ | 1 | 2.16330374133256 | 1.70597105956101 | +---+-------------------+---------------------+ | 2 | 2.02896192145109 | 1.61395085974027 | +---+-------------------+---------------------+ | 3 | 2.13040446092206 | 1.60944808047909 | +---+-------------------+---------------------+ | 4 | 2.05760815511591 | 1.60943791248579 | +---+-------------------+---------------------+ | 5 | 2.11206412847833 | 1.60943791243410 | +---+-------------------+---------------------+ | 6 | 2.07244694614844 | 1.60943791243410 | +---+-------------------+---------------------+