Math 640: EXPERIMENTAL MATHEMATICS Spring 2009 (Rutgers University) Webpage

http://www.math.rutgers.edu/~zeilberg/math640_09.html

Last Update: May 18, 2009.

Description

Experimental Mathematics used to be considered an oxymoron, but the future of mathematics is in that direction. In addition to learning the philosophy and methodology of this budding field, students will become computer-algebra wizards, and that should be very helpful in whatever mathematical specialty they'll decide to do research in.

We will first learn Maple, and how to program in it. This semester we will explore the fascinating field of combinatorial statistical physics, and critical phenomena. We will explore rigorous, semi-rigorous, and non-rigorous approaches, mostly using symbolic computation, but also numeric computations, using "Monte Carlo". But the actual content is not that important, it is mastering the methodology of computer-generated and computer-assisted research that is so crucial for your future.

There are no prerequisites, and no previous programming knowledge is assumed. Also, very little overlap with previous years, so people who have taken it in previous years are more than welcome to take it again, and become even greater wizards. The final projects for this class may lead to journal publications.


Added March 18, 2009: Pick a final project.

Diary and Homework

Programs done on Thurs., Jan. 22, 2009

Homework for Thurs., Jan. 22, class (due Jan. 26, 2009)

  1. Read and do all the examples, plus make up similar ones, in the first 30 pages of Frank Garvan's awesome Maple booklet.
  2. Do a first reading of the "textbook", just to get a feel for the subject.

Programs used and done on Monday, Jan. 26, 2009

Homework for Monday, Jan. 26, 2009, class (due Jan. 29, 2009)

  1. (for novices only) Read and do all the examples, plus make up similar ones, in pages 30-60 of Frank Garvan's awesome Maple booklet.
  2. (mandatory for experts, highly recommended for novices) Modify everything in Jan26.txt, so that the program B(n) becomes B(n,p), where now you toss a loaded coin with probability of winning a dollar being p (and losing a dollar being 1-p). Normalize the walk by first subtracting the expectation n(2p-1), and then dividing by the standard deviation sqrt(4np(1-p)). Of course you first must change fc(n) to lc(n,p) for a loaded coin. (Hint: if p=a/b, then use rand(1..b)() to generate an integer between 1 and a and decide that it is -1 if you get something between 1 and a and 1 if you get something between a+1 and b) [corrected Jan. 29, 2009, thanks to Dennis Hou].
  3. (mandatory for experts, highly recommended for novices) Write a program G(n,k), that inputs non-neg. integers n and k and outputs the set of sequences of {-1,1} of length n that add up to k and such that all the partial sums are non-negative. For example G(2,0); should yield {[1,-1]} . G(3,0); should yield {}, and G(4,0); should yield {[1,-1,1,-1],[1,1,-1,-1]} .

Programs done on Thursday, Jan. 29, 2009

Homework for Thurs., Jan. 29, 2009, class (due Feb. 2, 2009)

  1. (for novices only) Read and do all the examples, plus make up similar ones, in pages 60-90 of Frank Garvan's awesome Maple booklet. Hand-in a sample of 2 pages.
  2. (For everyone!) Using the methods of experimental mathematics, conjecture a closed-form expression for gnk(n,k) (note that if n+k is odd, then it is 0), the number of ways of starting at the origin, walking n steps, and winding-up at k, where at each step you can go forward or backward one unit. Once you conjectured it, prove it rigorously, by induction, using the "defining recurrence" of gnk(n,k) (plus the boundary conditions and initial conditions). Deduce, as a corollary, the closed-form expression (2n)!/(n!(n+1)!) for gnk(2n,0) (the famous Catalan numbers)
  3. ($5 to be divided among all solvers by Monday). Find a neat combinatorial proof (w/o looking it up!) of the Chung-Feller phenomenon that the number of ways of walking 2n steps from the origin back to the origin, with exactly 2k steps that lie to the left of the origin, is independent of k, and hence they are all equal to each other for k=0,1,2, ..., n, and hence each of them has binomial(2n,n)/(n+1) members, in particular proving once again that gnk(2n,0) is given by the Catalan numbers. [Added Feb. 2, 2009: Dennis Hou won this prize!].
  4. (for everyone) Modify the Gambler's ruin simulation GRe, in Jan29a.txt, to take another input, p, a rational number, where p is the chance of winning a dollar.
  5. (For everyone) In class we conjectured, using simulation, that the probability of a gambler (with a fair coin), who currently has i dollars, and continues to play until he or she lose everything or has n dollars, is i/n, and the expected duration of the game is i(n-i). Prove these facts rigorously and humanly. (Hint: set-up a system of linear equations for p(i), the prob. of winning with i dollars (and n fixed), and analogously for f(i), the duration of the game).

Programs done on Mon., Feb. 2, 2009

Homework for Mon., Feb. 2 2009, class (due Feb. 5, 2009)

  1. (for novices only) Read and do all the examples, plus make up similar ones, in pages 90-120 of Frank Garvan's awesome Maple booklet. Hand-in a sample of 2 pages.
  2. In class we empirically found out that the asymptotics for the Catalan numbers is 4n/n3/2/Pi1/2. Prove it rigorously, using Stirling's approximation for n!, or otherwise.
  3. Experiment with
    [seq(W(i,i,S), i=1..50)];
    and
    [seq(G(i,i,S), i=1..50)];
    for various choices of S. Use Zinn to conjecture asymptotics μn nθ with "nice" (rational) θ. Also experiment with the "probability of being good"
    [seq(G(i,i,S)/W(i,i,S), i=1..50)];
    and its asymptotics. Is θ always the same?
  4. Another way of thinking of W(m,n,S) is as the number of walks in the 2D square lattice from the origin to [m,n] using the step sizes from S. For G(m,n,S) it is the number of such walks that stay in m ≥ n . Write analogous programs for 3D, W3D(m,n,k,S) and G3D(m,n,k,S) (that stay in m ≥ n ≥ k. ). Check if the sequences
    [seq(W3D(i,i,i,S), i=1..20)];
    and
    [seq(G3D(i,i,i,S), i=1..20)];
    are in Sloane for S={1}, S={1,2}. Experiment with other S's and apply Zinn to them. Also for the sequence of probabilities of being good
    [seq(G3D(i,i,i,S)/W3D(i,i,i,S), i=1..30)];

Programs done on Thurs., Feb. 5, 2009

Homework for Thurs., Feb. 5 2009, class (due Feb. 9, 2009)

  1. Modify LD(p,x) to write a procedure LD2(p,x,y) that inputs a polynomial p in both x and y, with the meaning that the coeff. of xiyj is the probability of lending the pair [i,j], and outputs a random throw. (We would need this when we will study random walks in the plane).
  2. [Challenging!] Write a procedure GRP(p,x,N) that inputs a polynomial p in x (describing, as above, a certain [usually] loaded die), and a positive integer N, and outputs the list of length N-1 whose i-th entry is the prob. of winning in the GR game if you start at i.
    (Hint: Set up a system of linear equations and a set of variables and use the Maple solve command.)
  3. Assuming that you did the above challenging procedure, test it by comparing the first output of GGR for various i's and N's.
  4. [equally Challenging!] Write a procedure GRP(p,x,N) that inputs a polynomial p in x (describing, as above, a certain [usually] loaded die), and a positive integer N, and outputs the list of length N-1 whose i-th entry is the expected duration of the GR game if you start at i.
    (Hint: Set up a system of linear equations and a set of variables and use the Maple solve command.)

Programs done on Mon., Feb. 9, 2009

Homework for Mon., Feb. 9 2009, class (due Feb. 12, 2009)

  1. The corrected version of SymProbs has the line
    eq:={seq(subs(i=j,GS)=0,j=ldegree(P,x)+1..0)} union {seq(subs(i=j+N,GS)=1,j=0..degree(P,x)-1)}:
    which works. The previous version, in class, was
    eq:={seq(subs(i=j,GS)=0,j=ldegree(P,x)..0)} union {seq(subs(i=j+N,GS)=1,j=0..degree(P,x))}:
    Explain why it didn't work before, but it works now.
  2. Test the corrected procedure SymProbs(p,x,i,N), by comparing its output (with the given numeric N) to Probs(p,x,N), for the following polynomials p and pos. integers N
  3. ED(p,x,N), gives the expected durations of the game, no matter what. Modify it and write two procedures: EDlose(p,x,N) and EDwin(p,x,N) that computes the duration conditioned on the fact that you lost and won respectively. Check that their outputs add-up to the output of ED(p,x,n)
  4. (challenging!) Modify SymProbs(p,x,i,N) to write a program SymED(p,x,i,N) that outputs a formula for the expected duration of the game.
    Hint: if x=1 is a double root (the fair case) then you first find a particular solution that is quadratic in i, i.e.
    PS=ai^2+ bi +c
    plug-it into the recurrence,expand, and compare coefficients of powers of i, and equate it to 1, and solve for a, b,c. If x=1 is not a double root, then you try
    PS=ai+ b
    then the new general solution is
    PS+ WhatWeHadBefore

Programs done on Thu., Feb. 12, 2009

Special Homework for Thu. Feb. 12, 2009, class (due Feb. 14, 2009) [modified Feb. 14, 2009, 11:10am] (No extension!)


Added Feb. 14, 2009: A contest for the nicest picture, a chocolate box.
For c any complex mumber, define a transformation Mc(z)=z^2+c. A complex number c is (K,L)-good if iterating the map Mc K times, starting at z=0, i.e., doing z0=0 and zi+1=Mc(zi), for i=0, ..., K-1, in other words McK(0), has absolute value less than L.
  1. Write a procedure IsGood(c,K,L) that inputs a complex number c (in the form c=a+b*I, where a and b are floating-points) and outputs true or false according to whether c is (K,L)-good or not, respectively.
  2. Write a program MS(a,R,K,L) that inputs a number a, a small positive real number R (the resolution), and pos. integer K and L (K ≥ 10, L ≥ 3), and outputs a picture of all the (K,L) good points, with resolution R, i.e. all the points of the form A*R+B*R*I, with A and B integers such that A*R and B*R are both between -a and a. in the complex plane.
    Hints 1. to get the point corresponding to the complex number c, do [coeff(c,I,0),coeff(c,I,1)]. 2. Use plot(S,style=point, axes=none); to plot a set of points.
  3. Put the code in a file and call it vd. Then download the following input file invd, (keeling the name invd).
    Experiment with several choices of a,R,K,L, and pick the nicest picture to hand-out. The smaller the R, the better picture, but it will take much longer. Then, type
    maple < invd
    and look out for a file called vd.gif at the same directory. Print it out
(Added Feb. 16, 2009: the winners of the contest were Kellen Myers's animated Mandelbrot set (but he cheated, using Mathematica) and Emilie Hogan's pretty Mandelbrot set. They each get a chocolate box as a prize.)

Regular Homework for Thu. Feb. 12, 2009, class (due Feb. 16, 2009)

  1. (Modification of last time's challenging problem).
    $10 award to be divided among all correct solutions.
    Modify SymProbs(p,x,i,N) to write a program SymED(p,x,i,N) that outputs a formula for the expected duration of the game.
    Hint: if x=1 is a double root (the fair case) then you first find a particular solution that is quadratic in i, i.e.
    PS=ai^2
    plug-it into the recurrence,
    PS-add(coeff(P,x,j)*subs(i=i-j,PS),j=ldegree(P,x)..degree(P,x))=1
    expand, and compare coefficients and solve for a. If x=1 is not a double root, then you try, as above, but with
    PS=ai
    then the new general solution is
    PS+ WhatWeHadBefore
    and proceed as before. Now test your output by plugging-in specific values for N and i against the output of ED(p,x,N);
  2. The binary trees we talked about in class are a little weird in the sense that it is possible to have a left son but no right son, and vice versa. A complete binary tree is an ordered tree where every vertex has either no children, or exactly two children. Write a Maple program CBT(k) that inputs a pos. integer k, and outputs all complete binary trees with exactly k leaves. For example CBT(1)={[]}; CBT(2)={[[],[]]}; CBT(3)={[[],[[],[]]],[[[],[]],[]]}. Conjecture an exact formula for nops(CBT(k)).
  3. $2 dollars to be divided: prove your conjecture about nops(CBT(k)).
  4. Modify the drawing program DT to draw these kinds of trees. Note that now we don't have {} in the data structure, just []'s.

Programs done on Mon., Feb. 16, 2009

Homework for Mon. Feb. 16, 2009, class (due Feb. 19, 2009)

  1. Modify CBT(k) to write a program CTT(k) that inputs a pos. k and outputs all the complete ternary trees (ordered trees such that each vertex has either 0 children or 3 children) with exactly k leaves.
  2. Modify CBT(k) to write a program CQT(k) that inputs a pos. k and outputs all the complete quad trees (ordered trees such that each vertex has either 0 children or 4 children) with exactly k leaves.
  3. Modify cbt(k) to write a program ctt(k) that inputs a pos. k and outputs the number of complete ternary trees (ordered trees such that each vertex has either 0 children or 3 children) with exactly k leaves. Check whether [seq(ctt(k), k=1..20)] is in Sloane. Using Sloane, or otherwise, conjecture c closed form formula, and try and prove it.
  4. Modify cbt(k) to write a program cqt(k) that inputs a pos. k and outputs the number of complete quad trees (ordered trees such that each vertex has either 0 children or 4 children) with exactly k leaves. Check whether [seq(cqt(k), k=1..20)] is in Sloane. Using Sloane, or otherwise, conjecture c closed form formula, and try and prove it.
  5. Modify RCBT(p,k,Sad, Happy) to write a program RCTT(p,k,Sad,Happy) that does the analogous think but regarding complete ternary trees.
  6. Using human ingenuity, try to find a formula, in terms of p, for the probability of a random complete ternary family tree lasting for ever. Test your formula by running RCTT many times.

Programs done on Thur., Feb. 19, 2009

Homework for Thur. Feb. 19, 2009, class (due Feb. 23, 2009)

  1. At the beginning of today's class we showed how to derive everything "analytically" (or rather algebraically) about branching processes. For the complete binary tree, the probability generating function, F=F(p,z), (over all finite ("dead") trees) satisfies the equation
    F=(1-p)*z+p*z*F2 .
    By plugging-in z=1, you can get an equation for the probability of dying out. By doing subs(z=1,diff(F,z)) you get the expected size.
    Do the analogous things for complete ternary trees and complete quad trees, and verify that while the "boiling point" (transition prob.) are different (namely 1/2,1/3,1/4 resp.) the critial exponents of the quantities (obtained by setting p=pc+x, and looking at the lowest-order term in x), are all the same.
  2. (challenge, 5 dollars, I don't know the answer myself) Experiment with CER(n,eps,K) for fairly large (but not too large!) inputs of n and K. For eps not 0, it should be close to 1. When I did it for eps negative I got something way too small. Is it the program's fault? Were Erdos and Renyi wrong? Did Slade have a typo? Is it correct, but the n and K that we can do in real time are too small? I'd like to know.
  3. A triangle in a graph are three vertices {i,j,k} such that {{i,j}, {j,k}, {i,k}} are edges. Write a program NT(G,n) that inputs a graph (as a set of pairs) of {1, ..., n} and outputs the number of triangles.
  4. Can you find the "transition point", for p, for the number of triangles, when it goes from few and far between to lots and lots?
  5. A tree is a connected graph on n vertices without cycles (and it is easy to see that it must have n-1 edges). Write a program RG1(n,p) that does the same as RG(n,p) but quits as soon as it has n-1 edges, or returns FAIL. Write another program that runs RG1(n,p) many times (say K), and counts how many times you get a tree and how many times you don't. Can you guess a transition probability that goes from being unlikely to likely?

Programs done on Mon., Feb. 23, 2009

Homework for Mon. Feb. 23, 2009, class (due Feb. 26, 2009)

  1. (Given last time, but prize raised to 10 dollars, I don't know the answer myself) Experiment with CER(n,eps,K) for fairly large (but not too large!) inputs of n and K. For eps not 0, it should be close to 1. When I did it for eps negative I got something way too small. Is it the program's fault? Were Erdos and Renyi wrong? Did Slade have a typo? Is it correct, but the n and K that we can do in real time are too small? I'd like to know.
    [Added Feb. 26, 2009: Congratulations to Liyang Diao for figuring out the discrepancy. She won the $10 prize. See her writeup and read the Erdos-Renyi original 1961 article, kindly found by Liyang Diao.
    [Added March 2: Read Gordon Slade's Email clarifying the apparent inaccuracy (posted with his kind permission)]
  2. Write a program SAWd(n,d) that inputs two pos. integers, n and d, and outputs the set of self-avoiding walks in the d-dimensional (usual, hyper-cubic) lattice.
  3. Write a program SAWp(n) that inputs a positive integer n, and outputs the set of self-avoiding walks on n steps, that stay in the positive quadrant (x &ge 0, y ≥ 0).
  4. Write a program SAWab(n,a,b) that inputs a positive integer n, and non-neg. integers a and b, and outputs the set of self-avoiding walks on n steps, that start at the origin, that stay in the -a ≤ y &le b.

Programs done on Thurs., Feb. 26, 2009

Homework for Thur. Feb. 26, 2009, class (due March 2, 2009)

  1. Thanks to Michael Ratner, we know that Zinn (and also our own FCmt(L)) gives reasonable answers if you extract the even parts and odd parts separately. (to find the real μ, you have to take the square-root of the μ for the even (or odd) parts, of course, but theta is the same. Dig out from Sloane the sequences for the number of self-avoiding walks, in the plane, for the so-called triangular lattice, and for the so-called hexagonal (honey-comb) lattice, and see whether you get the same theta. Do you have to do the separation also here? Perhaps you need to do it mod 3? Experiment!
  2. The μ for the hexagonal lattice is conjectured by physicists to be sqrt(2+sqrt(2)). Check whether this is reasonable.
  3. Find the μ and θ for the number of self-avoiding walks in three dimensions, using Zinn and using FMCt(L) (compare notes).
  4. Let F(n,a,b,S) be the number of self-avoiding walks from the origin [0,0] to the point [a,b] never visiting the members of the set of point S. We are primarily interested in F(n,a,b,{}) but need the extra parameter S to facilitace an induction. Write a program with option remember to compute F(n,a,b,S).
    (Hint F(n,a,b,S)=F(n-1,a-1,b,S union {[a,b]})+F(n-1,a+1,b,S union {[a,b]})+F(n-1,a,b-1,S union {[a,b]})+F(n-1,a,b+1,S union {[a,b]}) ), also set up the initial condition and boundary conditions: if abs(a)+abs(b)>n then it is 0, if member([a,b],S) then it is 0, and if member([0,0],S) then it is 0. )
    Express s(n), the total number of self-avoiding walks in the square lattice in terms of F(n,a,b,{}), and write a program to compute s(n). Compare it with nops(SAW(n)); What takes longer to compute?

Programs done on March 2.,2009 (Virtual Class, Self-Study!)

Today, March 2, 2009, is a snow day, but I expect you to study on your own. Please read and understand carefully the new programs in the file Mar2.txt, that contains the following procedures

Homework for Mon. March 2, 2009, class (due March 5, 2009)

  1. It is conjectured that the expected distance-squared from the beginning (the origin), of a self-avoiding walk of length n, in the 2D square-lattice grows like Cn2 ν, where ν is one of the famous critical exponents, and C is some constant. Using the enhanced random-saw-generator sRSAW(n,K,C), with, say, K=5, C=2000 (or whatever you find often does not FAIL), write a program AvD(n,K,C,M) that runs sRSAW(n,K,C) M times (M failry big, say, at least a few hundreds), and takes the average of the distance-squared of the last-visited-point of all of them (that didn't FAIL). Do it for many n's (say 50 ≤ n ≤ 100). Try to estimate ν empirically. (Warning: sRSAW(n,K,C) is not guaranteed to generate a saw uniformly-at-random, so your answer may disagree with the one conjectured by physicists, but let's hope for the best.
  2. Modify sn(n) (and its subroutine snab(n,a,b), and its subsubroutine FabSn(a,b,S,n), to write procedures snR(n,c,d) (and subroutine snabR(n,a,b,c,d), and subsubroutine FabSnR(a,b,S,n,c,d), that counts the number of n-step self-avoiding walks on the 2D square-lattice, confined to the horizontal strip c ≤ y ≤ d, for any integers c and d satisfying c ≤ 0 ≤ d.
  3. Test your program by typing
    [seq(snR(n,0,1),n=0..14)];
    and look it up in Sloane. Who is the famous mathematician who scooped you?
  4. Load into a Maple session the Salvy/Zimmermann Maple package gfun, by typing
    with(gfun):
    and use guessgf to guess a generating function for the above mentioned sequence.
  5. Try to look-up in Sloane
    [seq(snR(n,c,d),n=0..9)];
    for various other c ≤ 0 ≤ d. Are they there yet?

Programs and Papers discussed on March 5,2009 (Guest Lecturer: Prof. Ilias Kotsireas)

Homework for Mon. March 2, 2009, class (due March 5, 2009)

  1. Experiment with several other cubic (and quadratic!) and even quartic polynomials, and get beautiful color pictures.
  2. Write a program F(m) to generate the combinatorial fractal discussed in the paper, using the following steps
  3. Implement the homotopy method for solving the cubic in Dr. Kotsireas' writeup.

Programs done on Mon. Mar. 9, 2009

Homework for Mon. March 9, 2009, class (due March 12, 2009)

  1. Experiment with gSAW(n,a,b,c,d,PAT) for various a,b,c,d, PAT. Find out if any of them are in Sloane. If they are long enough apply Zinn from GuessHolo2, and estimate μ and θ . Also try guessgf from gfun.
  2. Another interesting set of combinatorial objects, reminisicent of self-avoiding walks, are cube-free words on two letters, and square-free words on three letters. A word in a given alphabet contains a square if it has a stutter, i.e. a sequence of letters followed by an exact replica. For example: neither "stutter" nor "iloveloveyou" are squre-free, but "doron" is. Write a procedure SqFr(n) that inputs a pos. integer n and outputs the set of words (lists) of length n in the alphabet {0,1,2} that are square-free. Look up [seq(nops(SqFr(n)),n=1..10)] in Sloane. (Hint: first write a procedure called IsBad(w) that tests whether it has a "square" at the end, then build the sets recursively).
  3. Do the analogous thing for cube-free words on two letters {0,1}.

Programs done on Thur. Mar. 12, 2009

Homework for Thur. March 12, 2009, class (due March 23, 2009)

  1. Write a procedure IsGood(G,M,N) that inputs a graph (given as a set of edges) and returns true if and only if the connected component containing the vertex [0,0] has at least one vertex of the form [M,j] for some j.
  2. Write a procedure, MC(M,N,p,K), that generates K random subgraphs of [0,M]x[0,N] (using RP(M,N,p)) and counts the ratio of good ones.
  3. By taking K fairly big, experiment with p=1/10,2/10, ..., 9/10, and try to estimate the prob. θ(M,N,p) of a random graph being good (percolating)
  4. What we did so far, was bond percolation. Write analogous programs for site percolation where you choose to include, or exclude any vertex.
  5. Have a good Spring Break!

Added March 18, 2009: Pick a final project.

Programs done on Mon. Mar. 23, 2009

Homework for Mon. March 23, 2009, class (due March 26, 2009)

  1. Modify MC(a,b,p,x,K) to write a program EL(a,b,p,x,K) that empirically estimates the expected duration of
  2. Modify F(a,b,p,x) to write a program EL(a,b,p,x,K) that computes EXACTLY the expected duration of such a game.
  3. Using Zinn and Findrec of GuessHolo2, study the sequence 1/2-F(i,i,p,x), for p=x+x2, for p=x+x2+ ...+ x6, and for the two-dice Backgammon, with the Backgammon convention.
  4. (5 dollars to first solver): Using Maple, or even by hand, find a closed-form formula for the generating function
    sum(sum(F(a,b,p,x)*X^a*Y^b,a=0..infinity),b=0..infinity)
    where p=(x+x2)/2. Can you generalize it for for general dice? Check your answers against F(a,b,p,x).

Programs done on Mon. Mar. 26, 2009

Homework for Thurs. March 26, 2009, class (due March 30, 2009)

  1. Explain why this system gives the generating functions for N(a,b) and P(a,b)
    Eq:={N+(x+x^2)/2*P-1/(1-x)/(1-y)+x+x^2/2, P+(1/2)*(y+y^2)*N-1/(1-x)/(1-y)+1+y/2};
  2. Generalize FindGF(x,y) to write a procedure FindGFg(x,y,p,x) that automatically computes the generating function for the probabibilty of winning if it is your turn to move and you are a units away and your opponent is b units away and the prob. dist. of the die is given by the polynomial p of x. FindGFg(x,y,(x+x^2),x); should be the same as FindGF(x,y);
  3. Consider a Solitaire game, and let E1(a) be the expected number of moves until you get to the negative region, with an arbitrary die p(x). Write a Monte Carlo program that estimates the number of needed moves.
  4. As above, but a program that computes the EXAXT expectation
  5. Write a procedure FindGFE(x,p,x) that outputs the exact generating function (as a rational function of x)
    f(x):=Sum(E1(a)*x^a,a=0..infinity)

Programs done on Mon., Mar. 30, 2009

Homework for Monday, March 30, 2009, class (due April 2, 2009)

  1. Michael Ratner pointed out correctly that our definition of "IsGood" (percolating) does not correspond to the usual one, that only requires that there is any paths of 1's from the bottom to the top. Modify IsGood(M) to write a procedure IsOK(M) that returns true iff there is a paths os such 1's. (Hint: you can definte a procedure "NeighborOf" and iterate it).
  2. Write MCok(A,B,K,p) the analog of MC(A,B,K,p) for the more lenient definition of percolating.
  3. Write a procedure Histogram(A,B,K,r) that runs MC(A,B,K,p) for (K fairly big) and p=i/r where r is fairly big, for i=0...r, and outputs a histogram (discrete plot) where the height above p=i/r is MC(A,B,K,p). Experiment with various A,B, and r's .
  4. Write a procedure HistogramOK(A,B,K,r) that runs MCok(A,B,K,p) for (K fairly big) and p=i/r where r is fairly big, for i=0...r, and outputs a histogram (discrete plot) where the height above p=i/r is MC(A,B,K,p). Experiment with various A,B, and r's .
  5. Write a procedure TotalProb(A,B,p) that for (rather small A and B) finds the EXACT probability, expressed as a polynomial in p of a matrix being good.
  6. Write a procedure TotalProbOK(A,B,p) that for (rather small A and B) finds the EXACT probability, expressed as a polynomial in p of a matrix being OK. Plot it for various A and B's (not too big!)

Programs done on Thurs. April 2, 2009

Homework for Thurs.,April 2, 2009, class (due April 6, 2009)

  1. Rederive empirically (using Pade in the numapprox package) the exact expression for
    sum(nops(States(i))*x^i,i=1..infinity)
  2. Prove it rigorously, using computer and/or human ingenuity.

    [Added April 13, 2009: Congratulations to Dennis Hou for finding a beautiful human/computer proof.]

  3. Write a procedure IsComp(v1,v2) that outputs two states v1 and v2 (of the same length), i.e. vectors of {0,1,2} with the property that no 1 can be next to a 2, and outputs true if and only it is possible that v2 lies under v1 in a percolating matrix.
  4. Verify empirically Eqs. (2) (3) and (4) in Herb Wilf's paper "On the zeros of the Riesz' function in the Analytic Theory of Numbers

Programs done on Monday April 6, 2009

Homework for Mon.,April 6, 2009, class (due April 13, 2009)

  1. Write a program GFP(N,p,z) that inputs a pos. integer N and a varible z, and outputs the rational function
    Sum(PerPf(M,N,p)*z^M,M=1..infinity)
    (Hint: first write GPF1(M,s,z), for computing the generating funcions
    Sum(PerPf1(M,N,p,s)*z^M,M=1..infinity)
    for all states s in States(N) ( you need to set a system of linear equations mimicking PerPf1, but computing them all at once, since they need each other).
  2. Plot PerPf1(N,N,p) for N=2,3,4,... as far as your computer permits
  3. (10 dollars): Solve with Maple (if possible), American Mathematical Monthly Problem Number 11426 (April 2009, p. 365) to find
    (GAMMA(1/14)*GAMMA(9/14)*GAMMA(11/14))/
    (GAMMA(3/14)*GAMMA(5/14)*GAMMA(13/14))

    [Added April 13, 2009: Congratulations to Dennis Hou for winning the prize (and finding a Maple bug!). Read Dennis Hou's brilliant solution.]

The class on Thur. April 9, 2009 (Guest Lecturer: Andrew Baxter)

Many thanks to Andrew Baxter for subbing on Passover!

Homework for Thurs.,April 9, 2009, class (due April 13, 2009)

  1. Try plotting each of the 5 equations/functions listed at the bottom of the page at http://math.rutgers.edu/~baxter/StudySheets/Graphs.html
  2. Try to complete one of the Maple labs listed at http://math.rutgers.edu/~baxter/251-S08/index.html#Labs (check the gray boxes). Use any listed data you want, or make up your own.

Programs done on Monday April 13, 2009

We used some of the ideas in Persi Diaconis' fascinating paper .

Homework for Mon., April 13, 2009, class (due April 16, 2009)

  1. Using the list of English words (given as lists in lower-case letters) construct the 27 by 27 matrix whose rows and columns are indexed by the 26 letters in the Latin alphabet, as well as the space character, construct the matrix M (given as a list of lists) such that M[x][y] is the frequency of letter x followed by letter y in the list English. (Every letter x at the beginning of a word contibutes to M[Space,x] and every letter x at the end of a word contibutes to M[x,Space], while x followed y in a word contibutes to M[x,y].
  2. Write a program that inputs a text in English (in the 27-letter alphabet, ignoring punctuation and capitalization), viewed as one long word w=[w1,w2, ...] and a permutation, pi, of 27 lettters, and outputs its encrypted version [pi[w1],pi[w2], ...]
  3. (Challenge) Using the Metropolis algorithm as described in Persi Diaconis' fascinating article, write a program that decrypts a message encrypted by a simple substitution method. Experiment with "Mary had a little lamb her fleece was white as snow" and other sentences of your choice, by first encrypting them with a random permutation, and then decrypting them with Metropolis.

Programs done on Thurs., April 16, 2009

Homework for Thurs., April 16, 2009, class (due April 20, 2009)

  1. [Same as the challenge last time, but now it is a contest. $10 for the fastest program that would decipher a random encrypted message] Using the Metropolis algorithm as described in Persi Diaconis' fascinating article, write a program that decrypts a message encrypted by a simple substitution method. Experiment with "Mary had a little lamb her fleece was white as snow" and other sentences of your choice, by first encrypting them with a random permutation, and then decrypting them with Metropolis. Use the 27 by 27 matrix of the relative frequency of consecutive-letter-pairs (use S for space), (that you were supposed to find last time), where the neighbor of a permutation is all its 27*26/2 permutations obtained by transposing two entries ) with the "importance weight" of a permutation pi being (w is the encrypted word that you have to decrypt)

    mul(M[pi[w[i]]][pi[w[i+1]]],i=1..nops(w)-1):

    [Note added April 20, 2009: It turns out that the table is not faithful to real life. Here is Kellen Myers' better table. The contest is postponed to the last day of classes, May 4, 2009]
    [Note Added April 21, 2009: Here is Liang Diao's better table.

  2. Use the Metropolis algorithm to write a procedure RM(m,n,z,w,K) that inputs pos. integers m and n, pos. real numbers z and w, and a large pos. integer K, and outputs a "typical" {-1,1} m by n matrix with the parameters z and w, obtained by starting with a random m by n {-1,1} matrix, and applies the Jump K times.
  3. Write a program PlotM(M) that inputs a {-1,1} matrix and outputs a picture of M where -1 is a white square and 1 is a black square.
  4. Write a procedure MHCor(m,n,z,w,p1,q1,p2,q2) to empirically estimates the correlation of the random variables A[p1][q1] and A[p2][q2] (what Cor(m,n,z,w,p1,q1,p2,q2) does) but that will work for large m and n. Test it also for small values of m and n against the exact results of Cor(m,n,z,w,p1,q1,p2,q2) for various numerical choices of z and w.

Programs done on Mon., April 20, 2009

Homework for Mon., April 20, 2009, class (due April 23, 2009)

  1. Write a procedure that inputs a matrix A and decides whether it is "ordered" or not (it is ordered if the number of sign-changes is small). By using GenM(20,20,z,1,3000), try to find, empirically, the "critical point" in z, where it passes from disordered to ordered.
  2. Write a a procedure RN(A) that inputs a matrix A[i][j] of dimensions m by n, say (both divisible by 3) and outputs an m/3 by n/3 matrix such that B[i,j] is the "majority vote" of the nine entries {A[3*i+a,3*j+b], a=0,1,2, b=0,1,2}.
  3. Write a procedure RNG(k,z,w) that inputs k, uses GenM(3k,3k,z,w,3000) to generate a typical Ising matrix with parameters z,w and applies RN(A) k-1 times. Call (z,w) good if you get the all 1 or all -1 3 by 3 matrix, and call it bad, if you don't. Empirically, find the "phase transition" z=f(w), for any given w (from w=1 to w=3 with intervals of .1) where it moves from good to bad (or vice-versa).
  4. Write the analog of Ising1(z,w,t) for 2 by n matrices, getting a certain rational function in z,w,t.

Programs done on Thurs., April 23, 2009

Homework for Thurs., April 23, 2009, class (due April 27, 2009)

  1. [$10] Fix OnsagerN(m,nu) and/or LarsN(m,nu) so that they would give the same outputs for the same inputs. [Added April 27, 2009: Brian Nakamura fixed it, and won the prize, see below]
  2. Adapt Isingm(z,w,t,m) so that we are talking about a cylinder
  3. Adapt Ising(z,w,t,m) so that we are talking about a torus. [Hint: you can either use the present approach with beginning and ending column, or set-up an 2m by 2m matrix and find its eigenvalues.

Programs done on Mon., April 27, 2009

Brian Nakamura corrected OnsagerN (cylindrical convention), here is the the corrected version of Apr23.txt.

It was a nice day, so we used the 14 computers between our shoulders, as well as some input/output device called paper-and-pencil (or pen). We discussed how to write a webbook "120 KenKen puzzles", using team effort. The challenge is to have such a webbook ready to be posted in my homepage by Monday, May 4, 2009 (the last day of classes). The project leader is Andrew Baxter, in charge of putting everything together, and the html part. He has to boss three team leaders

Additional Homework for Mon., April 27, 2009, class (due April 30, 2009)

In addition to the class-project, please do the following
  1. Try to empirically guess the phase-transition in what we call ν by defining a function f[N](nu)=log(OnsagerN(N,nu))/N, and plotting f[N] for large N (say N=40 and N=100) for ν between 0 and 1 (using plot). Can you guess the location of the phase-transition (hint: it should for large N tend to log(1+sqrt(2))/2=0.44068... (the root of sinh(2 ν)=1, see Eq. (5.12))
  2. Extend LarsN(m,nu) to LarsNg(m,nu,w) that takes w different than 1. By defining a function g[N](nu,w)=log(LarsN(N,nu,w))/N, plot it for N=4 or N=5 and for various w, and try to see how the location of the "phase transition" (where it starts shooting up) depenes on w.

Programs done on Thu., April 30, 2009

Andrew Baxter kindly supplied the following Onsager plot clearly showing the phase transition at ν=log(1+sqrt(2))/2=0.4406867934... .

Homework for Thur. April 30, 2009, class (due May 4, 2009)

  1. Define a generalized directed animal with base k to be a trapezoid of 0's and 1's with the same rule as a directed animal, but with the top row consisting of all k ones (and no zeros!) . Write a program NuTk(n,k): that counts the number of such creatures with n ones.

    Example:

      1 1     1 1      1 1
     1 0 1   1 1 0    1 0 0
    
    are directed animals with base 2, and
      1 1
     1 0 0
    1 1 1 0
    
    would have been, but is not, since the third 1 at the bottom row has both its parents 0.
      1 1 1
     1 0 1 1
    0 1 1 0 1
    
    is a directed animal with base 3, etc.
  2. Write a program NuTall(n) that finds the total number of generalized animals in other words:
    NuTall:=proc(n): add(NuTk(n,k),k=1..n): end:
    Conjecture an explicit solution to NuTall(n) as an expression in n. [Added May 4, 2009: the answer is 3n-1 ].
  3. ($50 if you don't peek at the literature): prove your conjecture.

What we did on Monday, May 4, 2009

Today was the last day of classes, and we started with a deciphering competition. The students were supposed to use their Metropolis-based (simple substitution encryption) decription to decipher the following coded message. The winner was Kellen Myers. See his program (coming up soon, right now the link does not exist). Kellen won the book "Graph Theory" by Robin Wilson.

The second winner was Liyang Diao, who only finished a few seconds later. The answer was the first paragraph of this webpage. Liyang won a granola bar.

The second contest was to do as many Amer. Math. Montly problems, from the March, April, and May, 2009 issues, (that were handed out) using Maple, as possible, in the remaining 70 minutes. The winner was Eric Rowland who solved (and corrected!) Amer. Math. Monthly problem 11435. Eric won the book "A Walk through combinatorics" by Miklos Bona.

Many other students came close to solving several other problems. They are welcome to finish them up at home.

Homework for May 4, 2009, class

This was the last day of classes. Don't forget to Email me (at least a preliminary version of) your final project. Have a great summer!
Look at the students' Final Projects Note: these will be hopefully expanded and upgraded through the summer and beyond.
Added May 19, 2009: Look at the awsome collective class project Dr.Z.'s ExpMath Class's Calcudoku WebBook
Dr. Z.'s teaching page