<jupyter_start><jupyter_text>13.1 Apply Newton's method to the function $f(x) = (1-x)^2 + 100(1-x^2)^2=0$, using the initial guess of $2.5$.   Solution In order to use Newton's method, we must know the derivative of $f(x)$: $$f'(x) = 400 x^3 - 398x - 2.$$ $f(x)$ and $f'(x)$ are then defined using a $\texttt{lambda}$ function. The Newton's method solver, $\texttt{newton}$, found in Chapter 13 is then ran and the root approximated. The python file $\texttt{ch13.py}$ contains all of functions found in Chapter 13 and is imported for them to be used. It is also acceptable to paste the needed functions manually.<jupyter_code># Define f(x) and f'(x)
f = lambda x: (1-x)**2 + 100*(1-x**2)**2
f1 = lambda x: 400*x**3 - 398*x - 2

# Import functions from chapter 13
from ch13 import *

# Solve for root and print
root = newton(f,f1,2.5)
print('The root found is',root)<jupyter_output><empty_output><jupyter_text>Roots of the Bessel FunctionConsider the Bessel function of the first kind defined by $$ J_\alpha(x) = \sum_{m=0}^\infty \frac{(-1)^m}{m! \, \Gamma(m+\alpha+1)} {\left(\frac{x}{2}\right)}^{2m+\alpha}.  $$ - Write a Python code that prompts the user asks if they want to use bisection or Newton's method. Then the user enters an initial guess or initial range depending on the method selected. - Using the input find a root to $J_0$.  Each iteration print out to the user the value of $J_0(x)$ for the current guess and the change in the guess.- For testing, there is a root at $x \approx 2.4048$. Also, {\tt math.gamma(x)}, will give you $\Gamma(x)$.  Solution First, a function is defined that returns the value of $J_\alpha(x)$ to be used as an input for the root-finding methods. Note that $$\frac{d}{dx} J_\alpha(x) = -J_{\alpha+1}(x),$$ therefore we must also define a function for $-J_1(x)$ as the input for the derivative for Newton's method.$\br$ As the definition for both cases includes an infinite integral, it is a reasonable approximation to instead use a finite integral to a reasonably large value.$\br$<jupyter_code>import numpy as np
import math

# Define a function for J_0
def J_0(x, M = 100):
    """Order zero Bessel function of the first-kind
    evaluated at x
    Inputs:
        alpha:  value of alpha
        x:      point to evaluate Bessel function at
        M:      number of terms to include in sum
    Returns:
        J_0(x)
    """
    total = 0.0
    for m in range(M):
        total += (-1)**m/(math.factorial(m)*math.gamma(m+1))*(0.5*x)**(2*m)
    return total 

# Define a function for J_0prime
def J_0prime(x, M = 100):
    """Derivative of the order zero Bessel function
    of the first-kind evaluated at x
    Inputs:
        alpha:  value of alpha
        x:      point to evaluate Bessel function at
        M:      number of terms to include in sum
    Returns:
        J_0'(x)
    """
    total = 0.0
    for m in range(M):
        total += -(-1)**m/(math.factorial(m)*math.gamma(m+2))*(0.5*x)**(2*m+1)
    return total<jupyter_output><empty_output><jupyter_text>$\br$The $\texttt{bisection}$ function from Chapter 12 is copied and modified to print the results between iterations.$\br$<jupyter_code>def bisectionNew(f,a,b,epsilon=1.0e-6):
    """Find the root of the function f via bisection where the root lies within [a,b]
    Args:
        f: function to find root of
        a: left-side of interval
        b: right-side of interval
        epsilon: tolerance
    Returns:
        estimate of root
    """

    assert (b>a)
    assert (f(a)*f(b) < 0)
    delta = b - a
    print("We expect",int(np.ceil(np.log(delta/epsilon)/np.log(2))),"iterations")
    iterations = 0
    cOld = 0
    print('\n    f(guess)\t    Change in guess') # table header
    while (delta > epsilon):
        c = (a+b)*0.5
        change = cOld - c
        if (f(a)*f(c) < 0):
            b = c
        elif (f(b)*f(c) < 0):
            a=c
        else:
            return c
        delta = b-a
        iterations += 1
        if(iterations == 1):
            print(iterations,"  %.6e" % f(c)) # print initial guess
        else:
            print(iterations,"  %.6e" % f(c),"       %.6e" % change) # print guess w/ change
        cOld = c
    print("\nIt took",iterations,"iterations")
    return c #return midpoint of interval<jupyter_output><empty_output><jupyter_text>$\br$The $\texttt{newton}$ function from Chapter 13 is also copied and modified to print the results between iterations.$\br$<jupyter_code>def newtonNew(f,fprime,x0,epsilon=1.0e-6, LOUD=False):
    """Find the root of the function f via Newton-Raphson method
    Args:
        f: function to find root of
        fprime: derivative of f
        x0: initial guess
        epsilon: tolerance
        
    Returns:
        estimate of root
    """
    x = x0
    if (LOUD):
        print("x0 =",x0)
    iterations = 0
    print('\n    f(guess)\t  Change in guess') # table header
    while (np.fabs(f(x)) > epsilon):
        xOld = x
        if (LOUD):
            print("x_",iterations+1,"=",x,"-",f(x),"/",fprime(x),"=",x - f(x)/fprime(x))
        x = x - f(x)/fprime(x)
        iterations += 1
        change = x - xOld
        print(iterations,"  %.6e" % f(xOld),"      %.6e" % change) # print iteration
    print("\nIt took",iterations,"iterations")
    return x #return estimate of root<jupyter_output><empty_output><jupyter_text>$\br$A simple function, $\texttt{findRoot}$ is then defined that asks the user for the desired root-finding method. It then asks the user for the necessary inputs for each root-finding method and executes the necessary function.$\br$<jupyter_code>def findRoot():
    """Finds a root of the order zero Bessel function
    using either Newton's method or the method of
    bisection
    Args:
       None, inputs are done using input function 
    Returns:
        estimate of root
    """
    # Ask for method
    method = input('Enter a method (newton/bisection): ')
    
    # Newton's method
    if method == "newton":
        
        # Ask for initial guess
        print("\nUsing Newton's method")
        guess = float(input('Enter an initial guess: '))
        
        # Calculate the root
        root = newtonNew(J_0,J_0prime,guess)
    
    # Bisection method
    elif method == "bisection":
    
        # Ask for bounds
        print("Using bisecton method")
        a = float(input('\nEnter a bound to the left: '))
        b = float(input('Enter a bound to the right: '))
        
        # Calculate the root
        root = bisectionNew(J_0,a,b)
    
    # User you wrong
    else:
        
        print(method,"is not a proper method")
        
    return root<jupyter_output><empty_output><jupyter_text>$\br$First, we will test it with the method of bisection.$\br$<jupyter_code># Test with bisection
rootBi = findRoot()
print('\nThe result is',"%.6f" % rootBi)<jupyter_output><empty_output><jupyter_text>$\br$And second, the function is tested with Newton's method.$\br$<jupyter_code># Test with newton
rootNewton = findRoot()
print('\nThe result is',"%.7f" % rootNewton)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.12400906002265005, 'fraction_numerical': 0.018544733861834653, 'mean_word_length': 2.7963460505104782, 'pattern_counts': {'":': 2, '<': 32, '<?xml version=': 0, '>': 32, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '7', 'dir': 'jupyter-structured-clean-dedup', 'id': '271903', 'license': 'MIT', 'n_tokens_mistral': 2344, 'n_tokens_neox': 2258, 'n_words': 793, 'path': 'solution_chapter13.ipynb', 'repo_name': 'DrRyanMc/CompNucEng', 'size': '7064'}