Home
|
Course Outline
|
Lectures
|
Homework
|
Files
Here is the code which generates this applet:
// LatticeGas.java
import comphys.graphics.*;
import comphys.Easy;
import java.awt.*;
import java.awt.event.*;
public class LatticeGas extends Animation {
int[] lat = new int[10000];
int[][] nn = new int[10000][6];
int[] rule = new int[1024];
int[] nv = {3, 4, 5, 0, 1, 2};
int[] mask = {1, 2, 4, 8, 16, 32, 64, 128};
int lx = 17;
int ly = 20;
int bx = 1;
int by = 1;
int t; // time step
boolean flow; // inject particle at left edge
int injectDelay; // time steps to delay injection
void injectParticlesFromLeft () {
for (int j = 0; j < ly; j++) {
int n = j * lx;
lat[n] = 35;
}
}
void removeParticlesFromRight () {
for (int j = 0; j < ly; j++) {
int n = j * lx + lx - 1;
lat[n] &= 28;
}
}
double rho;
void computeDensity () {
rho = 0;
int sites = lx * ly;
for (int j = 0; j < ly; j++) {
for (int i = 0; i < lx; i++) {
int n = lat[j * lx + i];
if (n == 128) {
--sites;
} else {
for (int bit = 0; bit < 6; bit++) {
if ((n & 1) == 1)
rho += 1;
n >>= 1;
}
}
}
}
rho /= sites;
}
void initial () {
// begin with no particles
for (int n = 0; n < lx * ly; n++)
lat[n] = 0;
// fill block in center of lattice with 6 particles per site
for (int j = -by; j <= by; j++) {
for (int i = -bx; i <= bx; i++) {
int n = (j + ly / 2) * lx + i + lx / 2;
if (flow)
lat[n] = 128;
else
lat[n] = 63;
}
}
if (flow)
injectParticlesFromLeft();
computeDensity();
t = 0;
}
int[] latn = new int[10000];
void update () {
++t;
if (flow)
removeParticlesFromRight();
// bounce back boundary conditions v goes to -v
// move particles
for (int j = 0; j < ly; j++) {
for (int i = 0; i < lx; i++) {
int n = j * lx + i;
for (int dir = 0; dir < 6; dir++) {
if ((lat[nn[n][dir]] & mask[dir]) != 0) {
if ((lat[n] & mask[7]) != 0)
// reflection
latn[nn[n][dir]] |= mask[nv[dir]];
else
// particle moves from nearest neighbor
latn[n] |= mask[dir];
}
}
}
}
// collisions
for (int j = 0; j < ly; j++) {
for (int i = 0; i < lx; i++) {
int n = j * lx + i;
if ((lat[n] & mask[7]) == 0) {
lat[n] = rule[latn[n]];
latn[n] = 0;
}
}
}
if (flow) {
if (injectDelay == 0 || t % injectDelay == 0)
injectParticlesFromLeft();
}
}
void ruletable () {
// 6-bit saturated deterministic rule
for (int i = 0; i < 256; i++)
rule[i] = i;
rule[21] = 42;
rule[42] = 21;
rule[ 9] = 36;
rule[18] = 9;
rule[36] = 18;
rule[27] = 45;
rule[45] = 54;
rule[54] = 27;
rule[19] = 37;
rule[37] = 19;
rule[50] = 41;
rule[41] = 50;
rule[22] = 13;
rule[13] = 22;
rule[26] = 44;
rule[44] = 26;
rule[11] = 38;
rule[38] = 11;
rule[25] = 52;
rule[52] = 25;
}
void nntable () {
// nn[n][dir] gives the neighbor whose particle moving
// in direction dir will move to site n
for (int j = 0; j < ly; j += 2) {
for (int i = 0; i < lx; i++) {
int n = j * lx + i;
int ip = i + 1;
if (ip > lx - 1)
ip = 0;
int im = i - 1;
if (im < 0)
im = lx - 1;
int jp = j + 1;
if (jp > ly - 1)
jp = 0;
int jm = j - 1;
if (jm < 0)
jm = ly - 1;
nn[n][0] = j * lx + im;
nn[n][1] = jp * lx + im;
nn[n][2] = jp * lx + i;
nn[n][3] = j * lx + ip;
nn[n][4] = jm * lx + i;
nn[n][5] = jm * lx + im;
}
}
for (int j = 1; j < ly; j += 2) {
for (int i = 0; i < lx; i++) {
int n = j * lx + i;
int ip = i + 1;
if (ip > lx - 1)
ip = 0;
int im = i - 1;
if (im < 0)
im = lx - 1;
int jp = j + 1;
if (jp > ly - 1)
jp = 0;
int jm = j - 1;
if (jm < 0)
jm = ly - 1;
nn[n][0] = j * lx + im;
nn[n][1] = jp * lx + i;
nn[n][2] = jp * lx + ip;
nn[n][3] = j * lx + ip;
nn[n][4] = jm * lx + ip;
nn[n][5] = jm * lx + i;
}
}
}
boolean needToResize = true;
boolean makeBarrier;
class Lattice extends Plot implements MouseListener {
Lattice () {
addMouseListener(this);
}
double[] vx = new double[64];
double[] vy = new double[64];
double[] v = new double[64];
double[] radius = new double[64];
void createVelocityTables () {
// compute all possible combinations of velocities
double vmax = 0;
for (int i = 0; i < 64; i++) {
vx[i] = vy[i] = radius[i] = 0;
int n = i;
double theta = 0;
for (int j = 0; j < 6; j++) {
if ((n & 1) == 1) {
vx[i] += Math.cos(theta);
vy[i] += Math.sin(theta);
radius[i] += 1;
}
n >>= 1;
theta -= Math.PI / 3;
}
v[i] = Math.sqrt(vx[i] * vx[i] + vy[i] * vy[i]);
radius[i] = Math.sqrt(radius[i] / 12.0);
if (v[i] > vmax)
vmax = v[i];
}
// normalize them for plotting
for (int i = 0; i < 63; i++) {
if (v[i] > 0) {
vx[i] *= 1 / vmax;
vy[i] *= 1 / vmax * 2 / Math.sqrt(3);
}
}
}
double margin;
public void paint () {
clear();
double dx = 1;
double dy = Math.sqrt(3) / 2;
margin = (dy * (ly + 1) - dx * (lx + 1.5)) / 2;
if (margin < 0) { // fit horizontal
margin /= dy * dy;
setWindow(-1, lx + 0.5, margin - 1, ly - margin);
} else {
setWindow(-margin - 1, lx + 0.5 + margin, -1, ly);
}
if (needToResize) {
createVelocityTables();
needToResize = false;
}
for (int j = 0; j < ly; j++) {
for (int i = 0; i < lx; i++) {
int n = j * lx + i;
if (lat[n] == 128) {
double rx = 0.4;
double ry = rx * 2 / Math.sqrt(3);
setColor("red");
if (j % 2 == 0)
floodCircle(i - rx, i + rx, j - ry, j + ry);
else
floodCircle(i + 0.5 - rx, i + 0.5 + rx,
j - ry, j + ry);
} else if (lat[n] != 0) {
setColor("yellow");
double rx = radius[lat[n]];
double ry = rx * 2 / Math.sqrt(3);
if (j % 2 == 0)
floodCircle(i - rx, i + rx, j - ry, j + ry);
else
floodCircle(i + 0.5 - rx, i + 0.5 + rx,
j - ry, j + ry);
setColor("black");
if (v[lat[n]] != 0) {
if (j % 2 == 0)
plotLine(i, j,
i + vx[lat[n]], j + vy[lat[n]]);
else
plotLine(i + 0.5, j,
i + 0.5 + vx[lat[n]], j + vy[lat[n]]);
}
} else {
setColor("green");
if (j % 2 == 0)
plotPoint(i, j);
else
plotPoint(i + 0.5, j);
}
}
}
}
public void mouseClicked (MouseEvent me) { }
public void mouseEntered (MouseEvent me) { }
public void mouseExited (MouseEvent me) { }
public void mousePressed (MouseEvent me) { }
public void mouseReleased (MouseEvent me) {
int j = (int) Math.round(yWorld(me.getY()));
if (j < 0 || j >= ly)
return;
int i = 0;
if (j % 2 == 0)
i = (int) Math.round(xWorld(me.getX()));
else
i = (int) Math.round(xWorld(me.getX()) - 0.5);
if (i < 0 || i >= lx)
return;
int n = j * lx + i;
if (makeBarrier) {
if (lat[n] == 128)
lat[n] = 0;
else
lat[n] = 128;
} else {
if (lat[n] == 63)
lat[n] = 0;
else
lat[n] = 63;
}
repaint();
}
}
class Output extends Plot {
Output () {
setSize(400, 30);
}
public void paint () {
clear();
setWindow(0, 100, 0, 3);
setColor("blue");
boxArea(0, 100, 0, 3);
setColor("white");
plotString("t = " + t, 5, 1);
plotString("rho = " + Easy.format(rho, 4), 30, 1);
}
}
Lattice lattice;
Output output;
Reader lxReader, lyReader, bxReader, byReader;
Reader skipReader, flowReader, delayReader;
Checkbox makeBarrierCheckbox;
int skip;
public void init () {
initial();
nntable();
ruletable();
add(lattice = new Lattice());
add(output = new Output());
add(lxReader = new Reader("L_x = ", lx));
add(lyReader = new Reader("L_y = ", ly));
add(bxReader = new Reader("b_x = ", bx));
add(byReader = new Reader("b_y = ", by));
add(flowReader = new Reader("Flow ", flow));
add(delayReader = new Reader("Inject delay", injectDelay));
add(makeBarrierCheckbox = new Checkbox("Make barrier", makeBarrier));
makeBarrierCheckbox.addItemListener(new ItemListener() {
public void itemStateChanged(ItemEvent ie) {
makeBarrier = makeBarrierCheckbox.getState();
}
});
add(skipReader = new Reader("skip steps", skip));
addControlPanel();
}
public void step () {
for (int s = 0; s <= skip; s++)
update();
computeDensity();
lattice.repaint();
output.repaint();
}
public void reset () {
lx = lxReader.readInt();
ly = lyReader.readInt();
bx = bxReader.readInt();
by = byReader.readInt();
flow = flowReader.readBoolean();
injectDelay = delayReader.readInt();
skip = skipReader.readInt();
needToResize = true;
initial();
nntable();
ruletable();
lattice.repaint();
output.repaint();
}
public static void main (String[] args) {
LatticeGas latticeGas = new LatticeGas();
latticeGas.frame("Lattice Gas", 430, 680);
}
}
Questions or comments:
phygons@acsu.buffalo.edu