# AMath 583 Lab 5

## Where to find this notebook:

* This notebook will be in \$UWHPSC/labs/lab5/lab5a.ipynb if you have cloned the class repository.  You may need to "git pull" to update.
* Or you can bring it down to your computer (or SageMathCloud project) by typing this in a terminal shell:
  *   wget http://faculty.washington.edu/rjl/classes/am583s2014/lab5a.ipynb

## Announcements:

* Please sit at assigned tables.
* If you're having problems with notebooks in SMC, note that you can also use IPython from a terminal.

* Homework 2 will be posted by Thursday.  Part of the homework will be to read this paper: 
  "Best Practices for Scientific Computing, by G. Wilson, D. A. Aruliah", by C. T. Brown, et. al.
  that can be found at <http://www.plosbiology.org/article/info%3Adoi%2F10.1371%2Fjournal.pbio.1001745>
    


```
%pylab inline  
```

## SymPy

SymPy is a package of Python tools for symbolic calculations in algebra and calculus.  This is notebook has a very brief demo of some features, but you can learn much more about it from various on-line resources such as:

* <http://docs.sympy.org/latest/tutorial/intro.html#what-is-symbolic-computation>
* <http://nbviewer.ipython.org/github/jrjohansson/scientific-python-lectures/blob/master/Lecture-5-Sympy.ipynb>

You can do "from sympy import *" to get all of sympy in your namespace, but note that this will overwrite standard numpy variables such as pi and functions such as sin, cos, exp, sqrt, with sympy versions.  Since we may want to mix symbolic and numerical computing, we will be explicit about what's coming from SymPy:



```
import sympy as S  # shorthand
```


```
S.init_printing()   # so symbolic math is printed nicely
```

### Define a polynomial and factor it:


```
x = S.symbols('x')  # defines a symbol
f = x**3 - 3*x**2 - 18*x + 40   # a new symbolic expression
f                   # print it nicely
```


```
S.factor(f)
```

### We can also differentiate:


```
S.diff(f,x)  # differentiate f with respect to x
```

### A messier function and its derivative:


```
f = (x**3 * S.exp(5*x)*S.cos(S.pi*x)) / (1 + S.sqrt(x))
f
```


```
g = S.diff(f,x,n=3)
g
```

Note that if you "print g" it does not come out so pretty, but this might be more useful if you want to cut and paste this into a Fortran program, for example:


```
print g
```

### Evaluate this derivative at the point $x = 0.2$:


```
g2 = g.subs(x, 0.2)
g2
```

Note this substituted 0.2 for $x$, but the special symbol $\pi$ has not be converted to a number.  The SymPy function N converts special constants like $\pi$ and $e$ to numerical values:


```
S.N(g2)
```

You can specify how many digits to use when it does this substitution, e.g.


```
S.N(S.pi,n=50)
```

## Symbolic differentiation for a Newton iteration.

In the last two labs you implemented Newton's method for different functions and probably had to compute the derivatives by hand.  Here's a way to write and fvals function using SymPy:


```
def fvals(x, debug=False):
    from sympy import symbols,diff,sqrt
    # First specify f symbolically and differentiate using SymPy:
    xs = symbols('xs')
    fs = xs**2 - 2.
    fprimes = diff(fs, xs)
    
    # Now evaluate numerically at the value x passed in to this function:
    f = fs.subs(xs, x)
    fprime = fprimes.subs(xs, x)
    
    # The next lines are just for illustrating that this is working:
    if debug:
        print "fs = ",fs
        print "fprimes = ",fprimes
        print "x = ",x
        
    return f, fprime

# Try it out:
fv = fvals(3., debug=True)
print "fvals returns: ", fv
```

## Try this out with Newton's method:


```
def newton(fvals, x0, tol):
    xk = x0
    kmax = 30   # maximum number of iterations
    print "   k               xk                      f(xk)"
    for k in range(kmax):
        
        fxk, fpxk = fvals(xk)  # evaluate f(x) and f'(x)
        print "%4i   %22.15f  %22.15f"  % (k, xk, fxk)

        if abs(fxk) < tol:
            break  #leave the loop if tolerance satisfied
            
        xk = xk - fxk / fpxk  # update xk using Newton's method

    return xk
```


```
newton(fvals, 2, 1e-10)   
```

## Exercises:

* Use Newton's method to find a root of $f(x) = (x^2 - 2)\exp(-0.1 x^2) + 0.5$ from Lab 4.
* Let $f(x) = \sqrt{\cos(2\pi x) e^x}$.  Compute $f''(0.1)$, the second derivative evaluated at $x=0.1$.  This value is needed for the Lab 5 Quiz.

Remember you can get help in IPython by adding a ? to the end of an object and running the cell...


```
S.integrate?
```


```
S.diff?
```

819a1c4961d5cad46e6d2c6a9ea35edfd3cf4b04
ipynb
Jupyter Notebook
uwhpsc/labs/lab5/lab5a.ipynb
philipwangdk/HPC
e2937016821701adb80ece5bf65d43d1860640c0
MIT
uwhpsc/labs/lab5/lab5a.ipynb
philipwangdk/HPC
e2937016821701adb80ece5bf65d43d1860640c0
MIT
uwhpsc/labs/lab5/lab5a.ipynb
philipwangdk/HPC
e2937016821701adb80ece5bf65d43d1860640c0
MIT
Qwen/Qwen-72B
1. YES
2. YES