Differential equations and SageMath

The files below were on my teaching page when I was a college teacher. Since I retired, they disappeared. Samuel Lelièvre found an archived copy on the web, so I’m posting them here.

The files are licensed under the Attribution-ShareAlike Creative Commons license.

  1. Partial fractions handout, pdf
  2. Introduction to matrix determinants handout, pdf
  3. Impulse-response handout, pdf
  4. Introduction to ODEs, pdf
  5. Initial value problems, pdf
  6. Existence and uniqueness, pdf
  7. Euler’s method for numerically approximating solutions to DEs, pdf.
    Includes both 1st order DE case (with Euler and improved Euler) and higher order DE and systems of DEs cases, without improved Euler.
  8. Direction fields and isoclines, pdf
  9. 1st order ODEs, separable and linear cases, pdf
  10. A falling body problem in Newtonian mechanics, pdf
  11. A mixing problem, pdf
  12. Linear ODEs, I, pdf
  13. Linear ODEs, II, pdf
  14. Undetermined coefficients for non-homogeneous 2nd order constant coefficient ODEs, pdf
  15. Variation of parameters for non-homogeneous 2nd order constant coefficient ODEs, pdf.
  16. Annihilator method for non-homogeneous 2nd order constant coefficient ODEs, pdf.
    I found students preferred (the more-or-less equivalent) undetermined coefficient method, so didn’t put much effort into these notes.
  17. Springs, I, pdf
  18. Springs, II, pdf
  19. Springs, III, pdf
  20. LRC circuits, pdf
  21. Power series methods, I, pdf
  22. Power series methods, II, pdf
  23. Introduction to Laplace transform methods, I, pdf
  24. Introduction to Laplace transform methods, II, pdf
  25. Lanchester’s equations modeling the battle between two armies, pdf
  26. Row reduction/Gauss elimination method for systems of linear equations, pdf.
  27. Eigenvalue method for homogeneous constant coefficient 2×2 systems of 1st order ODEs, pdf.
  28. Variation of parameters for first order non-homogeneous linear constant coefficient systems of ODEs, pdf.
  29. Electrical networks using Laplace transforms, pdf
  30. Separation of variables and the Transport PDE, pdf
  31. Fourier series, pdf.
  32. one-dimensional heat equation using Fourier series, pdf.
  33. one-dimensional wave equation using Fourier series, pdf.
  34. one-dimensional Schroedinger’s wave equation for a “free particle in a box” using Fourier series, pdf.
  35. All these lectures collected as one pdf (216 pages).
    While licensed Attribution-ShareAlike CC, in the US this book is in the public domain, as it was written while I was a US federal government employee as part of my official duties. A warning – it has lots of typos. The latest version, written with Marshall Hampton, is a JHUP book, much more polished, available on amazon and the JHUP website. Google “Introduction to Differential Equations Using Sage”.

Course review: pdf

Love, War, and Zombies, pdf
This set of slides is of a lecture I would give if there was enough time towards the end of the semester

Integral Calculus and SageMath

Long ago, using LaTeX I assembled a book on Calculus II (integral calculus), based on notes of mine, Dale Hoffman (which was written in word), and William Stein. I ran out of energy to finish it and the source files mostly disappeared from my HD. Recently, Samuel Lelièvre found a copy of the pdf of this book on the internet (you can download it here). It’s licensed under a Creative Commons Attribution 3.0 United States License. You are free to print, use, mix or modify these materials as long as you credit the original authors.

Table of Contents

0 Preface

1 The Integral
1.1 Area
1.2 Some applications of area
1.2.1 Total accumulation as “area”
1.2.2 Problems
1.3 Sigma notation and Riemann sums
1.3.1 Sums of areas of rectangles
1.3.2 Area under a curve: Riemann sums
1.3.3 Two special Riemann sums: lower and upper sums
1.3.4 Problems
1.3.5 The trapezoidal rule
1.3.6 Simpson’s rule and Sage
1.3.7 Trapezoidal vs. Simpson: Which Method Is Best?
1.4 The definite integral
1.4.1 The Fundamental Theorem of Calculus
1.4.2 Problems
1.4.3 Properties of the definite integral
1.4.4 Problems
1.5 Areas, integrals, and antiderivatives
1.5.1 Integrals, Antiderivatives, and Applications
1.5.2 Indefinite Integrals and net change
1.5.3 Physical Intuition
1.5.4 Problems
1.6 Substitution and Symmetry
1.6.1 The Substitution Rule
1.6.2 Substitution and definite integrals
1.6.3 Symmetry
1.6.4 Problems

2 Applications
2.1 Applications of the integral to area
2.1.1 Using integration to determine areas
2.2 Computing Volumes of Surfaces of Revolution
2.2.1 Disc method
2.2.2 Shell method
2.2.3 Problems
2.3 Average Values
2.3.1 Problems
2.4 Moments and centers of mass
2.4.1 Point Masses
2.4.2 Center of mass of a region in the plane
2.4.3 x-bar For A Region
2.4.4 y-bar For a Region
2.4.5 Theorems of Pappus
2.5 Arc lengths
2.5.1 2-D Arc length
2.5.2 3-D Arc length

3 Polar coordinates and trigonometric integrals
3.1 Polar Coordinates
3.2 Areas in Polar Coordinates
3.3 Complex Numbers
3.3.1 Polar Form
3.4 Complex Exponentials and Trigonometric Identities
3.4.1 Trigonometry and Complex Exponentials
3.5 Integrals of Trigonometric Functions
3.5.1 Some Remarks on Using Complex-Valued Functions

4 Integration techniques
4.1 Trigonometric Substitutions
4.2 Integration by Parts
4.2.1 More General Uses of Integration By Parts
4.3 Factoring Polynomials
4.4 Partial Fractions
4.5 Integration of Rational Functions Using Partial Fractions
4.6 Improper Integrals
4.6.1 Convergence, Divergence, and Comparison

5 Sequences and Series
5.1 Sequences
5.2 Series
5.3 The Integral and Comparison Tests
5.3.1 Estimating the Sum of a Series
5.4 Tests for Convergence
5.4.1 The Comparison Test
5.4.2 Absolute and Conditional Convergence
5.4.3 The Ratio Test
5.4.4 The Root Test
5.5 Power Series
5.5.1 Shift the Origin
5.5.2 Convergence of Power Series
5.6 Taylor Series
5.7 Applications of Taylor Series
5.7.1 Estimation of Taylor Series

6 Some Differential Equations
6.1 Separable Equations
6.2 Logistic Equation

7 Appendix: Integral tables

Ring theory, via coding theory and cryptography

In these notes on ring theory, I tried to cover enough material to get a feeling for basic ring theory, via cyclic codes and ring-based cryptosystems such as NTRU. Here’s a list of the topics.

1 Introduction to rings
1.1 Definition of a ring
1.2 Integral domains and fields
1.3 Ring homomorphisms and ideals
1.4 Quotient rings
1.5 UFDs
1.6 Polynomial rings
1.6.1 Application: Shamir’s Secret Sharing Scheme
1.6.2 Application: NTRU
1.6.3 Application: Modified NTRU
1.6.4 Application to LFSRs

2 Structure of finite fields
2.1 Cyclic multiplicative group
2.2 Extension fields
2.3 Back to the LFSR

3 Error-correcting codes
3.1 The communication model
3.2 Basic definitions
3.3 Binary Hamming codes
3.4 Coset leaders and the covering radius
3.5 Reed-Solomon codes as polynomial codes
3.6 Cyclic codes as polynomial codes
3.6.1 Reed-Solomon codes as cyclic codes
3.6.2 Quadratic residue codes
3.6.3 BCH bound for cyclic codes
3.6.4 Decoding cyclic codes
3.6.5 Cyclic codes and LFSRs

4 Lattices
4.1 Basic definitions
4.2 The shortest vector problem
4.2.1 Application to a congruential PKC
4.3 LLL and a reduced lattice basis
4.4 Hermite normal form
4.5 NTRU as a lattice cryptosystem

Calculus on graphs

In these notes, I tried to cover enough material to get a feeling for “calculus on graphs”, with applications to sports rankings and the Friendship Theorem. Here’s a list of the topics.

1 . Introduction
2. Examples
3. Basic definitions
3.1 Diameter, radius, and all that
3.2 Treks, trails, paths
3.3 Maps between graphs
3.4 Colorings
3.5 Transitivity
4. Adjacency matrix
4.1 Definition
4.2 Basic results
4.3 The spectrum
4.4 Application to the Friendship Theorem
4.5 Eigenvector centrality and the Keener ranking
4.6 Strongly regular graphs
4.7  Orientation on a graph
5. Incidence matrix
5.1 The unsigned incidence matrix
5.2 The oriented case
5.3 Cycle space and cut space
6. Laplacian matrix
6.1 The Laplacian spectrum
7  Hodge decomposition for graphs
7.1 Abstract simplicial complexes
7.2 The Bjorner complex and the Riemann hypothesis
7.3 Homology groups
8. Comparison graphs
8.1 Comparison matrices
8.2 HodgeRank
8.3 HodgeRank example

Gray codes

This is based on work done about 20 years ago with a former student Jim McShea.

Gray codes were introduced by Bell Labs physicist Frank Gray in the 1950s. As introduced, a Gray code is an ordering of the n-tuples in GF(2)^n = \{0,1\}^n such that two successive terms differ in only one position. A Gray code can be regarded as a Hamiltonian path in the cube graph. For example:

[[0, 0, 0], [1, 0, 0], [1, 1, 0], [0, 1, 0], [0, 1, 1], [1, 1, 1], [1, 0, 1], [0, 0, 1]]

These can be generalized to n-tuples of integers (mod m) in the obvious way.

Gray codes have several applications:

  • solving puzzles such as the Tower of Hanoi and the Brain [G],
  • analog-digital-converters (goniometers) [S],
  • Hamiltonian circuits in hypercubes [Gil] and Cayley graphs of Coxeter groups [CSW],
  • capanology (the study of bell-ringing) [W],
  • continuous space-filling curves [Gi],
  • classification of Venn diagrams [R],
  • design of communication codes,
  • and more (see Wikipedia).
brain-puzzle

The Brain puzzle

Here's a SageMath/Python function for computing Gray codes.
def graycode(length,modulus):
    """
    Returns the n-tuple reverse Gray code mod m.


    EXAMPLES:
        sage: graycode(2,4)
        
        [[0, 0],
         [1, 0],
         [2, 0],
         [3, 0],
         [3, 1],
         [2, 1],
         [1, 1],
         [0, 1],
         [0, 2],
         [1, 2],
         [2, 2],
         [3, 2],
         [3, 3],
         [2, 3],
         [1, 3],
         [0, 3]]

    """
    n,m = length,modulus
    F = range(m)
    if n == 1:
        return [[i] for i in F]
    L = graycode(n-1, m)
    M = []
    for j in F:
        M = M+[ll+[j] for ll in L]
    k = len(M)
    Mr = [0]*m
    for i in range(m-1):
        i1 = i*int(k/m)
        i2 = (i+1)*int(k/m)
        Mr[i] = M[i1:i2]
    Mr[m-1] = M[(m-1)*int(k/m):]
    for i in range(m):
        if is_odd(i):
            Mr[i].reverse()
    M0 = []
    for i in range(m):
        M0 = M0+Mr[i]
    return M0


REFERENCES

[CSW] J. Conway, N. Sloane, and A. Wilks, “Gray codes and reflection groups”, Graphs and combinatorics 5(1989)315-325

[E] M. C. Er, “On generating the N-ary reflected Gray codes”, IEEE transactions on computers, 33(1984)739-741

[G] M. Gardner, “The binary Gray code”, in Knotted donuts and other mathematical entertainments, F. H. Freeman and Co., NY, 1986

[Gi] W. Gilbert, “A cube-filling Hilbert curve”, Math Intell 6 (1984)78

[Gil] E. Gilbert, “Gray codes and paths on the n-cube”, Bell System Technical Journal 37 (1958)815-826

[R] F. Ruskey, “A Survey of Venn Diagrams“, Elec. J. of Comb.(1997), and updated versions.

[S] Web page of T. Sillke

[W] A. White, “Ringing the cosets”, Amer. Math. Monthly 94(1987)721-746

Sports ranking methods, 4

This is the fourth of a series of expository posts on matrix-theoretic sports ranking methods. This post discusses the Elo rating.

This system was originally developed by Arpad Elo (Elo (1903-1992) was a physics professor at Marquette University in Milwaukee and a chess master, eight-time winner of the Wisconsin State Chess Championships.) Originally, it was developed for rating chess players in the 1950s and 1960s. Now it is used for table tennis, basketball, and other sports.

We use the following version of his rating system.

As above, assume all the $n$ teams play each other (ties allowed)
and let r_i denote the rating of Team i, i=1,2,\dots,n.

Let A=(A_{ij}) denote an $n\times n$ matrix of score results:

A_{ij}= \left\{ \begin{array}{rr} -1,& {\rm if\ team\ } i {\rm \ lost\ to\ team\ } j,\\ +1,& {\rm if\ team\ } i {\rm\ beat\ team\ } j,\\ 0, & {\rm if}\ i=j. \end{array} \right.

Let S_{ij}=(A_{ij}+1)/2.

As in the previous post, the matrix A associated to the example of the Patriot league is the adjacency matrix of a diagraph.

  1. Initialize all the ratings to be 100: {\bf r}=(r_1,\dots,r_n) = (100,\dots,100).
  2. After Team i plays Team j, update their rating using the formula

    r_i = r_i+K(S_{ij}-mu_{ij}),

    where K=10 and

    \mu_{ij} = (1+e^{-(r_i-r_j)/400})^{-1}.

In the example of the Patriot league, the ratings vector is

{\bf r}=(85.124, 104.79, 104.88, 85.032, 94.876, 124.53).

This gives the ranking

Lafayette < Army < Lehigh < Bucknell < Holy Cross < Navy.

This gives a prediction failure rate of 13.3\%.

Some SageMath code for this:

def elo_rating(A):
    """
    A is a signed adjacency matrix for a directed graph.

    Returns elo ratings of the vertices of Gamma = Graph(A) 
        
    EXAMPLES:
        sage: A = matrix(QQ,[
        [0 , -1 , 1  , -1 , -1 , -1 ],
        [1,   0 ,  -1,  1,  1,   -1  ],
        [-1 , 1 ,  0 ,  1 , 1  , -1  ],
        [1 , -1 , -1,  0 ,  -1 , -1  ],
        [1 , - 1 , - 1 , 1 , 0 , - 1  ],
        [1 ,  1  ,  1  , 1  , 1  , 0 ]
        ])
        sage: elo_rating(A)
        (85.124, 104.79, 104.88, 85.032, 94.876, 124.53)

    """
    n = len(A.rows())
    RR = RealField(prec=20)
    V = RR^n
    K = 10
    r0 = 100 # initial rating
    r = n*[r0]
    for i in range(n):
        for j in range(n):
            if ij and A[i][j]==1:
                S = 1
            elif ij and A[i][j]==-1:
                S = 0
            else:
                S = 1/2
            mu = 1/(1+e^(-(r[i]-r[j])/400))
            r[i] = r[i]+K*(S-mu)
    return V(r)

Sports ranking methods, 1

This is the first of a series of expository posts on matrix-theoretic sports ranking methods. This post, which owes much to discussions with TS Michael, discusses Massey’s method.

Massey’s method, currently in use by the NCAA (for football, where teams typically play each other once), was developed by Kenneth P. Massey
while an undergraduate math major in the late 1990s. We present a possible variation of Massey’s method adapted to baseball, where teams typically play each other multiple times.

There are exactly 15 pairing between these teams. These pairs are sorted lexicographically, as follows:

(1,2),(1,3),(1,4), …, (5,6).

In other words, sorted as

Army vs Bucknell, Army vs Holy Cross, Army vs Lafayette, …, Lehigh vs Navy.

The cumulative results of the 2016 regular season are given in the table below. We count only the games played in the Patriot league, but not including the Patriot league post-season tournament (see eg, the Patriot League site for details). In the table, the total score (since the teams play multiple games against each other) of the team in the vertical column on the left is listed first. In other words, ”a – b” in row $i$ and column $j$ means the total runs scored by team i against team j is a, and the total runs allowed by team i against team j is b. Here, we order the six teams as above (team 1 is Army (USMI at Westpoint), team 2 is Bucknell, and so on). For instance if X played Y and the scores were 10-0, 0-1, 0-1, 0-1, 0-1, 0-1, then the table would read 10-5 in the position of row X and column Y.

X\Y Army Bucknell Holy Cross Lafayette Lehigh Navy
Army x 14-16 14-13 14-24 10-12 8-19
Bucknell 16-14 x 27-30 18-16 23-20 10-22
Holy Cross 13-14 30-27 x 19-15 17-13 9-16
Lafayette 24-14 16-18 15-19 x 12-23 17-39
Lehigh 12-10 20-23 13-17 23-12 x 12-18
Navy 19-8 22-10 16-9 39-17 18-12 x
sm261_baseball-ranking-application_teams-digraph

Win-loss digraph of the Patriot league mens baseball from 2015

In this ordering, we record their (sum total) win-loss record (a 1 for a win, -1 for a loss) in a 15\times 6 matrix:

M = \left(\begin{array}{cccccc} -1 & 1 & 0 & 0 & 0 & 0 \\ 1 & 0 & -1 & 0 & 0 & 0 \\ -1 & 0 & 0 & 1 & 0 & 0 \\ -1 & 0 & 0 & 0 & 1 & 0 \\ -1 & 0 & 0 & 0 & 0 & 1 \\ 0 & -1 & 1 & 0 & 0 & 0 \\ 0 & 1 & 0 & -1 & 0 & 0 \\ 0 & 1 & 0 & 0 & -1 & 0 \\ 0 & -1 & 0 & 0 & 0 & 1 \\ 0 & 0 & 1 & -1 & 0 & 0 \\ 0 & 0 & 1 & 0 & -1 & 0 \\ 0 & 0 & -1 & 0 & 0 & 1 \\ 0 & 0 & 0 & -1 & 1 & 0 \\ 0 & 0 & 0 & -1 & 0 & 1 \\ 0 & 0 & 0 & 0 & -1 & 1 \end{array}\right).

We also record their total losses in a column vector:

{\bf b}= \left(\begin{array}{c} 2 \\ 1 \\ 10 \\ 2 \\ 11 \\ 3 \\ 2 \\ 3 \\ 14 \\ 4 \\ 14 \\ 10 \\ 11 \\ 22 \\ 6 \\ \end{array}\right).

The Massey ranking of these teams is a vector {\bf r} which best fits the equation

M{\bf r}={\bf b}.

While the corresponding linear system is over-determined, we can look for a best (in the least squares sense) approximate solution using the orthogonal projection formula

P_V = B(B^tB)^{-1}B^t,

valid for matrices B with linearly independent columns. Unfortunately, in this case B=M does not have linearly independent columns, so the formula doesn’t apply.

Massey’s clever idea is to solve

M^tM{\bf r}=M^t{\bf b}

by row-reduction and determine the rankings from the parameterized form of the solution. To this end, we compute

M^tM= \left(\begin{array}{rrrrrr} 5 & -1 & -1 & -1 & -1 & -1 \\ -1 & 5 & -1 & -1 & -1 & -1 \\ -1 & -1 & 5 & -1 & -1 & -1 \\ -1 & -1 & -1 & 5 & -1 & -1 \\ -1 & -1 & -1 & -1 & 5 & -1 \\ -1 & -1 & -1 & -1 & -1 & 5 \end{array}\right)

and

M^t{\bf b}= \left(\begin{array}{r} -24 \\ -10 \\ 10 \\ -29 \\ -10 \\ 63 \\ \end{array}\right).

Then we compute the rref of

A= (M^tM,M^t{\bf b}) = \left(\begin{array}{rrrrrr|r} 5 & -1 & -1 & -1 & -1 & -1 & -24 \\ -1 & 5 & -1 & -1 & -1 & -1 & -10 \\ -1 & -1 & 5 & -1 & -1 & -1 & 10 \\ -1 & -1 & -1 & 5 & -1 & -1 & -29 \\ -1 & -1 & -1 & -1 & 5 & -1 & -10 \\ -1 & -1 & -1 & -1 & -1 & 5 & 63 \end{array}\right),

which is

rref(M^tM,M^t{\bf b})= \left(\begin{array}{rrrrrr|r} 1 & 0 & 0 & 0 & 0 & -1 & -\frac{87}{6} \\ 0 & 1 & 0 & 0 & 0 & -1 & -\frac{73}{6} \\ 0 & 0 & 1 & 0 & 0 & -1 & -\frac{53}{6} \\ 0 & 0 & 0 & 1 & 0 & -1 & -\frac{92}{3} \\ 0 & 0 & 0 & 0 & 1 & -1 & -\frac{73}{6} \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 \end{array}\right).

If {\bf r}=(r_1,r_2,r_3,r_4,r_5,r_6) denotes the rankings of Army, Bucknell, Holy Cross, Lafayette, Lehigh, Navy, in that order, then

r_1=r_6-\frac{87}{6},\ \ r_2=r_6-\frac{73}{6},\ \ r_3=r_6-\frac{53}{6},\ \ r_4=r_6-\frac{92}{6},\ \ r_5=r_6-\frac{73}{6}.

Therefore

Lafayette < Army = Bucknell = Lehigh < Holy Cross < Navy.

If we use this ranking to predict win/losses over the season, it would fail to correctly predict Army vs Holy Cross (Army won), Bucknell vs Lehigh, and Lafayette vs Army. This gives a prediction failure rate of 20\%.