Examples » Walker: Integrating the gamma SDE

Control file

title "Example problem"

walker

  term  35.0    # Max time
  dt    0.01    # Time step size
  npar  100000  # Number of particles
  ttyi  100     # TTY output interval

  rngs
    mkl_r250 seed 1 end
  end

  gamma
    depvar l
    init zero
    coeff const
    ncomp 2
    # k = bS/kappa, 1/theta = b(1-S)/kappa
    # <Y> = S/(1-S), <y^2> = kappa/b * <Y>/(1-S)
    b     1.5            2.5 end
    kappa 1.0            1.0 end
    S     0.666666666666 0.8 end
    rng mkl_r250
  end

  statistics
    <l1l1> <l2l2> <l1l2>
  end

  pdfs
    interval          100
    filetype          txt
    policy            overwrite
    centering         node
    format            scientific
    precision         4
    f( L1 : 2.0e-1 )
    g( L2 : 2.0e-1 )
  end
end

Example run on 8 CPUs

./charmrun +p8 Main/walker -v -c ../../tmp/gamma.q

Results

The left figure shows the time evolution of the means estimated from the numerical simulation together with those of the invariant distributions. The right figure shows the time evolution of the variances and those of the invariant. Both plots indicate convergence to the correct statistically stationary state.

Gnuplot commands to reproduce the above plots:

plot "stat.txt" u 2:3 w l t "<L1>", "stat.txt" u 2:4 w l t "<L2>", 2.0 lt 1, 4.0 lt 2
plot "stat.txt" u 2:5 w l t "<l1l1>", "stat.txt" u 2:7 w l t "<l2l2>", 4.0 lt 1, 8.0 lt 2

The left figure shows the 2 estimated PDFs at the final step of the simulation together with the corresponding invariants. The right figure shows the time evolution of the estimated covariance indicating no correlations at all times corresponding to the statistically independent equations integrated.

Gnuplot commands to reproduce the above plots:

plot "pdf_f.txt" t "k=1.0, theta=2.0", "pdf_g.txt" t "k=2.0, theta=2.0", x**(1.0-1.0)*exp(-x/2.0)/gamma(1.0)/2.0**1.0 lt 1, x**(2.0-1.0)*exp(-x/2.0)/gamma(2.0)/2.0**2.0 lt 2
plot "stat.txt" u 2:6 w l t "<l1l2>"