Difference between revisions of "CompSciWeek11"

From Predictive Chemistry
Jump to: navigation, search
(Created page with "= Class 1 = * Tensor Manipulations - distance distributions in random point set exercise ** [http://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.spatial.distance.pdi…")
 
 
(3 intermediate revisions by the same user not shown)
Line 1: Line 1:
  +
= Reading Assignment =
  +
* [http://apps.nrbook.com/empanel/index.html Numerical Recipes] Pages 37-40
  +
* [http://calnewport.com/blog/2012/10/26/mastering-linear-algebra-in-10-days-astounding-experiments-in-ultra-learning/ Feynman Learning Technique]
  +
** Try and identify the most challenging missing concept and explain it to someone else!
  +
* Optional Additional reference: Basis and Singular Value lectures from [http://codingthematrix.com/ codingthematrix.com]
  +
 
= Class 1 =
 
= Class 1 =
 
* Tensor Manipulations - distance distributions in random point set exercise
 
* Tensor Manipulations - distance distributions in random point set exercise
 
** [http://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.spatial.distance.pdist.html pdist]
 
** [http://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.spatial.distance.pdist.html pdist]
 
** [http://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.stats.histogram.html histogram]
 
** [http://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.stats.histogram.html histogram]
* Rodriguez' Rotation Formula - the preferred method when you ''have'' to work with an angle.
+
* Rodrigues' Rotation Formula - the preferred method when you ''have'' to work with an angle.
 
* Matrix Multiplication
 
* Matrix Multiplication
 
** Basic linear algebra operations
 
** Basic linear algebra operations
 
** Timing Strassen (see also Winograd)
 
** Timing Strassen (see also Winograd)
  +
  +
== Distribution Function Code ==
  +
  +
This illustrates dramatic simplifications that can be obtained by using tensors in numpy.
  +
  +
<source lang="python">
  +
# Radial distribution function calculator for un-scaled distributions.
  +
# (g(r) is get_rdf() * L**3/N and goes to 1).
  +
from numpy import *
  +
  +
# Get the unnormalized rdf in an isotropic cubic box of side length L
  +
# return unit is particles per L's length unit^3
  +
# Note that grid methods should be used for very large N.
  +
def get_rdf(x, L, M):
  +
D2 = x - x[:,newaxis]
  +
D2 -= L*floor(D2/L+0.5) # wrap into box
  +
r2 = sqrt(sum(D2*D2, -1))
  +
s = arange(M+1)*0.5*L/M
  +
h, _ = histogram(r2, s)
  +
h[0] -= len(x) # remove self-distances
  +
# Divide by volume (times an extra (N-1), since we counted
  +
# N*(N-1) distances, but the density scales as N)
  +
h = h.astype(float)
  +
h /= (len(x)-1) * 4*pi/3.0*(s[1:]**3 - s[:-1]**3)
  +
return h
  +
  +
# This should give a uniform RDF of 1000 pt / L**3 (= 1 here).
  +
def test():
  +
N = 1000
  +
L = 10.0
  +
x = random.random((N,3))*L
  +
print get_rdf(x, L, 20)
  +
</source>
  +
  +
Rodrigues' formula code. This is for completeness, so you can generate 3D rotations given in axis-angle notation. It is definitely not for memorization.
  +
  +
<source lang="python">
  +
# Build rotation matrix to rotate about an arbitrary vector using the right-hand rule.
  +
# Uses Rodrigues' rotation formula (in the quaternion representation).
  +
def build_rotation(u, theta):
  +
n = u
  +
m = sum(n*n) # check that n is normalized
  +
if(m < 1.0-1.0e-10 or m > 1.0+1.0e-10):
  +
n /= sqrt(m)
  +
  +
trans = zeros((3,3), float)
  +
  +
s = sin(theta)
  +
c = cos(theta)
  +
trans[0,0] = c + n[0]*n[0]*(1.0-c)
  +
trans[0,1] = n[0]*n[1]*(1.0-c) - n[2]*s
  +
trans[0,2] = n[1]*s + n[0]*n[2]*(1.0-c)
  +
  +
trans[1,0] = n[2]*s + n[0]*n[1]*(1.0-c)
  +
trans[1,1] = c + n[1]*n[1]*(1.0-c)
  +
trans[1,2] = -n[0]*s + n[1]*n[2]*(1.0-c)
  +
  +
trans[2,0] = -n[1]*s + n[0]*n[2]*(1.0-c)
  +
trans[2,1] = n[0]*s + n[1]*n[2]*(1.0-c)
  +
trans[2,2] = c + n[2]*n[2]*(1.0-c)
  +
  +
return trans
  +
</source>
   
 
= Class 2 =
 
= Class 2 =

Latest revision as of 18:12, 5 November 2014

Reading Assignment

Class 1

  • Tensor Manipulations - distance distributions in random point set exercise
  • Rodrigues' Rotation Formula - the preferred method when you have to work with an angle.
  • Matrix Multiplication
    • Basic linear algebra operations
    • Timing Strassen (see also Winograd)

Distribution Function Code

This illustrates dramatic simplifications that can be obtained by using tensors in numpy.

<source lang="python">

  1. Radial distribution function calculator for un-scaled distributions.
  2. (g(r) is get_rdf() * L**3/N and goes to 1).

from numpy import *

  1. Get the unnormalized rdf in an isotropic cubic box of side length L
  2. return unit is particles per L's length unit^3
  3. Note that grid methods should be used for very large N.

def get_rdf(x, L, M):

   D2 = x - x[:,newaxis]
   D2 -= L*floor(D2/L+0.5) # wrap into box
   r2 = sqrt(sum(D2*D2, -1))
   s = arange(M+1)*0.5*L/M
   h, _ = histogram(r2, s)
   h[0] -= len(x) # remove self-distances
   # Divide by volume (times an extra (N-1), since we counted
   #                   N*(N-1) distances, but the density scales as N)
   h = h.astype(float)
   h /= (len(x)-1) * 4*pi/3.0*(s[1:]**3 - s[:-1]**3)
   return h
  1. This should give a uniform RDF of 1000 pt / L**3 (= 1 here).

def test():

   N = 1000
   L = 10.0
   x = random.random((N,3))*L
   print get_rdf(x, L, 20)

</source>

Rodrigues' formula code. This is for completeness, so you can generate 3D rotations given in axis-angle notation. It is definitely not for memorization.

<source lang="python">

  1. Build rotation matrix to rotate about an arbitrary vector using the right-hand rule.
  2. Uses Rodrigues' rotation formula (in the quaternion representation).

def build_rotation(u, theta):

       n = u
       m = sum(n*n) # check that n is normalized
       if(m < 1.0-1.0e-10 or m > 1.0+1.0e-10):
               n /= sqrt(m)
       trans = zeros((3,3), float)
       s = sin(theta)
       c = cos(theta)
       trans[0,0] = c + n[0]*n[0]*(1.0-c)
       trans[0,1] = n[0]*n[1]*(1.0-c) - n[2]*s
       trans[0,2] = n[1]*s + n[0]*n[2]*(1.0-c)
       trans[1,0] = n[2]*s + n[0]*n[1]*(1.0-c)
       trans[1,1] = c + n[1]*n[1]*(1.0-c)
       trans[1,2] = -n[0]*s + n[1]*n[2]*(1.0-c)
       trans[2,0] = -n[1]*s + n[0]*n[2]*(1.0-c)
       trans[2,1] = n[0]*s + n[1]*n[2]*(1.0-c)
       trans[2,2] = c + n[2]*n[2]*(1.0-c)
       return trans

</source>

Class 2

  • Implicit solutions to linear algebraic equations
    • Least-squares fitting
    • Repetition for point sets - project out and re-fit
    • Review projection and orthogonality (see also Gram-Schmidt)
    • Solution by factorization + forward / reverse substitution
  • Minimization the hard way - nonlinear problems Optimize
    • Transcendental problem, 5th order polynomials, etc.
    • Geodesic paths via numerical optimization