Hybrid Path Integral Monte Carlo: Multiple-Grain Simulation

   Physics 
  This is a simulation of the Efetov-AES phase functional for granular metals.  
  Here there are just 4 identical grains arranged in a straight line,
  with periodic boundary conditions, labelled by X=0,1,2,3.
  Thus we have a 1+1D system.

  The kinetic term in the action couples "neighbouring beads on the same chain"
  (phase of a single grain at nearest-neighbour TAU slices).

  The dissipative term couples 4 beads:  2 on one chain and 2 on another chain.
  This interaction is long-ranged in TAU, but nearest-neighbour in X (grains).

  The simulation uses Hybrid Path Integral Monte Carlo,
  with Verlet time evolution.
  The Metropolis rejection step has been left out
  (and hence the algorithm does not obey detailed balance exactly,
  especially if ds is large!)

  
   Display 
  The upper panel displays the current phase configurations PHI(TAU,X)
  on the vertical axis running from 0 to 2*pi,
  versus imaginary time tau
  on the horizontal axis running from 0 to BETA.
  Phases on the four grains are plotted in black, red, green, blue resp.

  The lower panel currently displays nothing, but in future it can be used
  for displaying the imaginary-time electromagnetic response
  K(TAU) = dj/dA, which is related to the conductivity SIGMA(OMEGA).

  The status bar displays the MC step number.
  Note that each MC step consists of an MD trajectory
  (and a Metropolis rejection step, which is currently not implemented).
  The MD velocities are randomised at the beginning of every MD trajectory;
  this effect may be noticeable.

  
   Usage 
  The number of slices in imaginary time, JMAX, is fixed at 64.
  The following parameters may be changed by modifying the values in the boxes
  AND THEN PRESSING ENTER:

  m = mass = grain capacitance = C/e^2 = 1/E_c
  eta = damping = dissipation = inter-grain tunnelling conductance = Efetov's g
  beta = inverse temperature = 1/T
  ds = MD step size

  You probably want to set ds to be quite large (about 0.09),
  but not so large that the time evolution becomes unstable. 
  If eta is zero, then ds=sqrt(m/eps)/3 is a good choice.

  When BETA is small (high temperature), 
  quantum fluctuations of the path are small,
  and each phase behaves like a classical particle.

  When BETA is large (low temperature), 
  quantum fluctuations are large, 
  and non-zero winding numbers may start to appear spontaneously.

  An interesting regime is when BETA is large (50) but ETA is also large (10).
  Then, the phases experience a large drift with TAU, 
  but phases on adjacent grains at the same TAU are locked together.

*/ import java.awt.*; import java.awt.event.*; import javax.swing.*; public class MG extends JApplet implements Runnable, ActionListener, ComponentListener { //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% // Variable declarations final double pi=Math.PI,twopi=2*pi; final int iparams=20; final int xmax=4; // 4 grains long time_display=-0x7FFFFFFF,interval_display=10; long time_display_2=-0x7FFFFFFF,interval_display_2=1000; volatile Thread thread1; boolean running = true; JPanel jPanel0,jPanelL,jPanelR,jPanel1,jPanel2; JLabel[] jLabelArray; JTextField[] jTextFieldArray; MyPicture myPicture1,myPicture2; int g_Xmin,g_Xmax,g_Ymin,g_Ymax; double g_xmin,g_xmax,g_ymin,g_ymax; int jmax = 64; int cfgs; double tt[],uu[][],vv[][],ff[][],cc[],ss[],CC[],CCsum[],CCmean[],aa[]; double dd1[]; double m,kappa,eta,beta,eps,ds; //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% // Mathematical functions public double linmap (double x, double x1, double x2, double y1, double y2) { return y1 + (y2-y1)/(x2-x1)*(x-x1); } public double randomGaussian () { double x,y,u,t,xx,yy; do { x = Math.random()*2-1; y = Math.random()*2-1; u = x*x + y*y; } while (u>1); t = Math.sqrt(-2*Math.log(u)/u); xx = t*x; yy = t*y; return xx; } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% // Given displacements uu, compute forces ff // f[j] = m/eps*.. void computeForce() { int i,j,k; double dum1,dum2,dum3; int x,x2; for (x=0; xtime_display) { time_display = time_now + interval_display; for (x=0; xtime_display_2) { time_display_2 = time_now + interval_display_2; updatePicture2(); myPicture2.repaint(); } } computeForce(); for (x=0; x0) if (Math.abs(Y2-Y1)*2 < (g_Ymin-g_Ymax)) g.drawLine (X1, Y1, X2, Y2); X1 = X2; Y1 = Y2; } } } void updatePicture2() { Graphics g; Image image1; if ((image1=myPicture1.image1)==null) return; if ((g=image1.getGraphics())==null) return; g_Xmin=0; g_Xmax=myPicture2.getWidth(); g_Ymax=0; g_Ymin=myPicture2.getHeight(); g_xmin=0.0; g_xmax=beta; g_ymin=0.0; g_ymax=1.0; g.setColor (Color.white); clear (g); /* setColor (Color.blue); plotXY (g, jmax, tt, CCmean); */ } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% public void actionPerformed(ActionEvent e) { //if (e != null) { // && e.getActionCommand().equals ("p")) { readParamsFromTextFields(); } void readParamsFromTextFields() { double dum1; // Read primary parameters m = Float.parseFloat (jTextFieldArray[0].getText()); eta = Float.parseFloat (jTextFieldArray[1].getText()); beta = Float.parseFloat (jTextFieldArray[2].getText()); ds = Float.parseFloat (jTextFieldArray[3].getText()); // Recalculate secondary parameters and tables eps = beta/jmax; //ds = Math.sqrt(eps/m)/3; // /3 int j; for (j=0; j0) if (Math.abs(Y2-Y1)*2 < (g_Ymin-g_Ymax)) g.drawLine (X1, Y1, X2, Y2); X1 = X2; Y1 = Y2; } } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% // The MyPicture class contains an image. It is like a label. public class MyPicture extends JComponent { Image image1; protected void paintComponent(Graphics g) { g.drawImage (image1, 0, 0, null); } } }