• Generate a random number x according to the following distributions:

(a)                    2x exp(-x2) ,    x > 0
       p(x) =  {
                         0 ,    x < 0
 
  Solution:
 

        F(x) = int _{-infinity to x} p(x) dx

                 = int _{0 to x} 2x exp{-x² } dx
                 = int _{0 to x} exp{-x² } dx²
                 = 1 - exp(-x²)

          x = F^-1(xi) =  [ - ln( 1 - xi ) ]^1/2 = [ - ln xi ]^1/2
 

 

(b)                    exp(-x)/[ 1 + exp(-x) ] ,   x > 0
       p(x) = {
                         0 ,    x < 0
 
  Solution:
 

       F(x) = int _{-infinity to x} p(x) dx

                 = int _{0 to x} exp(-x) / [ 1 + exp(-x) ]  dx
                 = - int _{0 to x} 1/[ 1 + exp(-x) ] d exp(-x)
                 =  ln 2 - ln [1 + exp(-x) ]

          x = F^-1(xi) =  - ln( 2e^-xi - 1) 
 
 

(c)                     x exp(-x)   ,    x > 0
       p(x) = {
                          0  ,       x < 0

Solution:
 

       F(x) = int _{-infinity to x} p(x) dx

                 = int _{0 to x}  x exp{-x } dx
                 = - int _{0 to x} x d exp{-x }
                 = - x exp{-x } +  int _{0 to x}  exp{-x } dx
                 = 1 - ( x + 1 )exp(-x)

          xi = 1 - ( x + 1 ) exp(-x)

           This function can not be solved analytically. It needs to
            solved numerically. Write a C-code to do this.
 
 
 

• Write down the Monte Carlo estimator and the error for the following integrals:

(a)   int _{0 to infinity}  2x exp(-x2)
 

Solution:
 

Monte Carlo estimator:
 
 G = int   g(X) p(X) dX

G_N =  sum _{i=1 to N} g(Xi) / N

Monte Carlo error:

sigma² ~

 N/(N - 1){1/N sum _{i=1 to N} g²(Xi) -
                   (1/N sum _{i=1 to N} g(Xi))²}
 

g(x) = 2x

G_N =  sum _{i=1 to N} 2 Xi / N

sigma² ~
 N/(N - 1){1/N sum _{i=1 to N} 4 Xi² -
                   (1/N sum _{i=1 to N} 2Xi)²}
 
 

 
(b)   int _{0 to infinity} exp(-x)/[ 1 + exp(-x) ]
 
 
Solution:

g(x) = 1/ [ 1 + exp(-x) ]

G_N =  sum _{i=1 to N}  [ 1 + exp(-Xi) ]^-1 / N
 

sigma² ~
 N/(N - 1){1/N sum _{i=1 to N} [ 1+ exp(-Xi) ]^-2 -
                   (1/N sum _{i=1 to N} [ 1 + exp(-Xi) ]^-1)²}

 

(c)   int _{0 to infinity} x exp(-x)
 
 

Solution:
 

Monte Carlo estimator:
 
 g(x) = x

G_N =  sum _{i=1 to N}  Xi / N

sigma² ~
 N/(N - 1){1/N sum _{i=1 to N}  Xi² -
                   (1/N sum _{i=1 to N} Xi)²}

• The transition probability W(i to j)  between 3 states is:

           W(2 to 1) = 1/6, W(2 to 2) = 1/2,
           W(2 to 3) = 1/3

           W(3 to 1)=0, W(3 to 2)=2/3, W(3 to 3) = 1/3
 
        (a)  Show that the solution for the probabilities is:

               P1 = 1/10, P2 = 6/10, P3=3/10
 
        (b)  Use P0=(1,0,0) as the initial distribution,
               how fast the system converges to the equilibrium
               probability distribution?
 
 
          Solution:

          (a)

                       0       1       0

           W = 1/6     1/2     1/3

                      0     2/3     1/3
 
 
             P = ( P1, P2, P3)

 
             P = P W

             P1 =  0  + 1/6P2 + 0
             P2 = P1 + 1/2 P2 + 2/3 P3
             P3 = 0 + 1/3 P2 + 1/3 P3

              P1 + P2 + P3 = 1

              which gives:  P1 = 1/10, P2 = 6/10, P3=3/10
 

             (b)

              P0 = (1, 0, 0)
 
              P1 = P0 W = (0, 1, 0)

               P2 = P1 W = (1/6, 1/2, 1/6)

               P3 = P2 W = (1/12, 19/36, 2/9) = (0.08, 0.52, 0.22)

               P4 = P3 W = ....

• Consider a situation where we have only three possible states 1, 2, and 3 (e.g. discrete quantum states). The transition probability is given as

W(1 to 2) = 1, W(2 to 1) = W(2 to 3) = W(3 to 2) = W(3 to 3) = 1/2, and W(1 to 3) = W(3 to 1) = W(1 to 1) = W(2 to 2) = 0.
 

        (b) Calculate the equilibrium probability distribution Pi,
              where i = 1,2,,3.

        (c) Does the equilibrium distribution satisfy detailed balance
        condition?
 
 

Solution:

(a)

           0        1        0
 
 W =   1/2     0        1/2
 
            0       1/2     1/2
 
 

           1/2      0        1/2

W² =   0       3/4       1/4

            1/4    1/4       1/2
 
 

             0       3/4      1/4

W³ =    3/8    1/8      1/2

             1/8    1/2       3/8
 

 

W⁴  = nonzero matrix.

Thus the Markov chain is ergodic.
 

(b)

Equilibrium distribution P=(P1, P2, P3) satisfies the equation:

P = P W

which gives:

P1 =  0 + 1/2 P2 + 0
P2 = P1 + 0 + 1/2 P3
P3 = 0 + 1/2 P2 + 1/2 P3

In addition there is the condition:

P1+P2+P3 = 1

From these equations P can be solved.
 
 
 

(c)

Detailed balance condition:
 
           P(X') W(X' to X) = P(X) W(X to X')
 
 Check whether P1, P2, P3 satisfies:
 

 Pi W(i, j) = Pj W(j, i)
 

If yes, then P satisfies the detailed balance condition.