/* This C program uses the theory from Rosenthal (1995b) applied to the frog example to find a bound on the total variation distance after 200 iterations */ #define NODES 100 /* the size of the state space */ #define M 4 /* P has edges between nodes within M of each other */ #define ITERS 200 /* the number of jumps before we take an approximate sample from the stationary distribution */ #define BETA 2 /* the target distribution is proportional to e^(-((distance to node zero)^BETA)) */ #define HEXP 2 /* the exponent in the calculation of h below */ #include #include #include int distance(int); int main() { int min, i, j, x, y, disti, distj, distx, disty; double P[NODES][NODES], pi[NODES]; double bestbound=1, Eh0=0, Eh1=0, alph=100, epsilon, A=0, denom=0, acceptalpha, bound, temp; /* the following two loops calculate the stationary distribution of the chain */ for (i=0; i pi[i]) acceptalpha = 1; else { disti = distance(i); distj = distance(j); /* to avoid dividing zero by zero; this makes P slightly different than actually specified, but the difference is minor, and the same thing had to be done in the eigenvalue calculation */ if ((pi[i]disti)) acceptalpha = 0; else if ((pi[i]= NODES-M) { P[i][j] = acceptalpha/(2*M+1); P[i][i] = P[i][i] + (1-acceptalpha)/(2*M+1); } else if ((M+1 <= abs(i-j)) && (abs(i-j) < NODES-M)) P[i][j] = 0; } } /* finds expected value of h after 1 iteration and alph using h = 1 + dist(x,0)^HEXP + dist(y,0)^HEXP and taking the low value of temp since inequality holds for all x, y covering all pairs of x and y isn't necessary with this particular h since x will equal y, but it's nice to have the code in case we want to change h */ for (x=1; x