/** applet No. 1127 * * template - periodic molecular dynamics 3D * - Morse potential * - multi thread MorseMD :: display - asynchronous * * Created by Ikeuchi Mitsuru on February 24 2007. * Copyright (c) 2007 Ikeuchi Mitsuru. All rights reserved. * * ver 0.0.1 2007.02.24 created * ver 0.0.2 2007.06.17 improved code * ver 0.0.3 2007.07.08 change Particles class * */ import java.awt.*; import java.awt.event.*; import java.applet.*; public class MorseMtPMD3D extends Applet implements MouseListener, MouseMotionListener, ItemListener, ActionListener, AdjustmentListener, Runnable { // ------------------------------------ preset field ----------- Thread th = null; // for run()-paint() thread Particles ma = new Particles(); MolecularDynamics3D md = null; DisplayMD3D dsp = new DisplayMD3D(); // for event Choice ch_mat1, ch_mat2, ch_tempControl, ch_view; Button bt_reset, bt_startStop, bt_step; Scrollbar sc_temp, sc_scale; // for off-paint buffer Dimension dim; Image imgOff; Graphics gOff; int dgX,dgY, dgXb,dgYb; int sleepTime = 100; int dispMode = 0; int thCount = 0; // ------------------------------ applet main thread ----------- public void init() { resize(630,460); setBackground(Color.white); dim = getSize(); imgOff = createImage(dim.width,dim.height); gOff = imgOff.getGraphics(); ch_mat1 = new Choice(); for (int i=0; i="+(int)(md.meanPressure()/1.0e4+0.5)/100.0+"MPa",630/6*5+0,110); gOff.drawString("nkT="+(int)((Nmt/(xMax*yMax*zMax)*1.38e-23*tmp)/1.0e4+0.5)/100.0+"MPa",630/6*4+0,110); gOff.setColor(Color.black); gOff.drawString("regMax="+md.regMax()+"",630/6*5+10,130); gOff.drawString("secMax="+md.secMax()+"",630/6*5+10,150); gOff.setColor(Color.black); gOff.drawString("thread MD = "+md.getCount()+" ",20,440); gOff.drawString("thread disp = "+thCount+" ",220,440); } // ------------------------------------ plot methods ----------- private void ballPlot(int gbx, int gby, double size, int imode) { int i,j,k,jj,gx,gy,g2x,g2y,gr; double sz,sc,tm; setPos(); rotate(); zSort(); drawWaku(gbx,gby,0); sc = dispScale*viewScale; sz = size*sc; for (i=0; i0.99) tm = 0.99; gOff.setColor(Color.getHSBColor( (float)md.getColor(j),0.8f,(float)(0.8*tm+0.2))); gOff.fillOval(gx-gr/2,gy-gr/2, gr, gr); } } drawWaku(gbx,gby,1); } private double distance(int i,int j) { double x,y,z; x = px[i]-px[j]; y = py[i]-py[j]; z = pz[i]-pz[j]; return Math.sqrt(x*x+y*y+z*z); } private void velocityPlot(int gbx, int gby, double mag) { int i,j,gx,gy,g2x,g2y; double sc,m, p2x,p2y,p2z,rotp2x,rotp2y,rotp2z; setPos(); rotate(); zSort(); drawWaku(gbx,gby,0); sc = dispScale*viewScale; for (i=0; i=xMax/2.0) { gOff.setColor(Color.getHSBColor((float)md.getColor(j),0.8f,(float)(0.2*pz[j]/zDepth+0.5))); } else { gOff.setColor(Color.getHSBColor((float)md.getColor(j),0.4f,(float)(0.2*pz[j]/zDepth+0.5))); } gOff.drawOval(gx-gr/2,gy-gr/2, gr, gr); } drawWaku(gbx,gby,1); } // plot utility private void drawWaku(int gbx, int gby, int imode) { int i,gx,gy,g2x,g2y,farPoint; double sc,tmp; sc = dispScale*viewScale; farPoint = findMin(); for (i=0; i<12; i++) { gx = (int)(sc*pwkx[boxp[i][0]])+gbx; gy = (int)(sc*pwky[boxp[i][0]])+gby; g2x = (int)(sc*pwkx[boxp[i][1]])+gbx; g2y = (int)(sc*pwky[boxp[i][1]])+gby; if (imode==0) { if (boxp[i][0]==farPoint || boxp[i][1]==farPoint) { gOff.setColor(Color.lightGray); gOff.drawLine(gx,gy, g2x, g2y); } } else if (imode==1) { if (boxp[i][0]!=farPoint && boxp[i][1]!=farPoint) { gOff.setColor(Color.gray); if (md.getTempMode()==2) { tmp = md.temperature(); if (tmp>md.getContTemp()) { gOff.setColor(Color.black); } else { gOff.setColor(Color.red); } } gOff.drawLine(gx,gy, g2x, g2y); } } } } private int findMin() { int i,im; double m; im = 0; m = pwkz[im]; for (i=0; i<8; i++) { if (pwkz[i]=srtz[1][p]) p = i; k = qSortPartition(i,j,srtz[1][p]); qSort(i,k-1); qSort(k,j); } } private int qSortPartition(int i,int j, int x) { int l,r,w; l = i; r = j; while (l<=r) { while (l<=j && srtz[1][l]=i && srtz[1][r]>=x) r--; if (l>r) break; w = srtz[0][l]; srtz[0][l] = srtz[0][r]; srtz[0][r] = w; w = srtz[1][l]; srtz[1][l] = srtz[1][r]; srtz[1][r] = w; l++; r--; } return ( l ); } // rotate private void rotate() { int i; for (i=0; i0) { gOff.fillRect(gx+i*2,gy+hoganHight-iy,2,iy); } } drawHogan(gx, gy, hoganWidth, hoganHight, xDiv, yDiv); gOff.setColor(Color.black); gOff.drawString("0",gx-5,gy+hoganHight+20); gOff.drawString("1 (nm)",gx+hoganWidth-20,gy+hoganHight+20); gOff.drawString("radial distribution",gx+40,gy+hoganHight+40); } private void drawHogan(int gx, int gy, int hoganWidth, int hoganHight, int xDiv, int yDiv) { int i; gOff.setColor(Color.gray); for (i=0; i0) { gOff.drawLine(gx+i,gy+200-iy,gx+i,gy+200); } } // Maxwell if (mat<21) { setMaxwell( md.getContTemp(), mat ); gOff.setColor(Color.getHSBColor(0.01f,0.8f,0.6f)); for (i=0; i<200-1; i++) { iy = (int)(ascale*ffvMaxwell[i]); iy2 = (int)(ascale*ffvMaxwell[i+1]); if (iy>0) { gOff.drawLine(gx+i,gy+200-iy,gx+i,gy+200-iy2); } } // MAT1 gOff.setColor(Color.getHSBColor(md.getColorOf(mat),0.8f,0.6f)); for (i=0; i<200-1; i++) { iy = (int)(ascale*meanvdf[1][i]); iy2 = (int)(ascale*meanvdf[1][i+1]); if (iy>0) { gOff.drawLine(gx+i,gy+200-iy,gx+i,gy+200-iy2); } } // MAT2 //gOff.setColor(Color.getHSBColor(md.getColorOfKind(MAT2),0.8f,0.6f)); //for (i=0; i<200-1; i++) { // iy = (int)(ascale*meanvdf[2][i]); // iy2 = (int)(ascale*meanvdf[2][i+1]); // if (iy>0) { // gOff.drawLine(gx+i,gy+200-iy,gx+i,gy+200-iy2); // } //} } // grid gOff.setColor(Color.gray); for (i=0; i<200; i+=50) { gOff.drawLine(gx+i,gy, gx+i, gy+200); } for (i=0; i<200; i+=100) { gOff.drawLine(gx,gy+200-i, gx+200, gy+200-i); } gOff.setColor(Color.black); gOff.drawRect(gx,gy,200,200); gOff.setColor(Color.black); gOff.drawString("0",gx-5,gy+hoganHight+20); gOff.drawString("2000 m/s",gx+hoganWidth-30,gy+hoganHight+20); gOff.drawString("velocity distribution",gx+30,gy+hoganHight+40); gOff.setColor(Color.lightGray); gOff.drawString("count",gx+120,gy+20); if (mat<21) { gOff.setColor(Color.getHSBColor(md.getColorOf(mat),0.8f,0.6f)); gOff.drawString("mean distr.",gx+120,gy+40); //gOff.setColor(Color.getHSBColor(md.getColorOf(mat2),0.8f,0.6f)); //gOff.drawString("mean distr.",gx+120,gy+60); gOff.setColor(Color.getHSBColor(0.01f,0.8f,0.6f)); gOff.drawString("Maxwell(Tc)",gx+120,gy+80); } } private void setVelocityDistribution(int MAT1) { int i,iv,knd; double vv; for (i=0; i<220; i++) { vdf[0][i]=0; vdf[1][i]=0; vdf[2][i]=0; } for (i=0; i=0.8*rCutoff && r=0.8*rCutoff && r=600 && ir<1000) { rt = 0.5+0.5*Math.cos(3.1415926536*(((double)ir)/1000.0-0.6)/0.4); } else { rt = 0.0; } return( rt ); } // ----------------------------------- class methods ----------- public double force(double r, int ki, int kj) { int i; double a; i = (int)(r/hh); a = r - i*hh; return ( ((hh-a)*forceTable[ki][kj][i] + a*forceTable[ki][kj][i+1])/hh ); } // ------------------------------------- I/O methods ----------- public static int getKindMax() { return kindMax; } public static String getMatStr(int knd) { return matStr[knd]; } public static double getMassOf(int knd) { return ion[knd][0]; } public static double getDiaOf(int knd) { return 2.0*ion[knd][5]; } public static float getColorOf(int knd) { return (float)ion[knd][6]; } public static double getrCutoff() { return rCutoff; } public static double getdt() { return dt; } } */ // ============== molecular dynamics - Morse potential =========== class MolecularDynamics3D extends Thread { private Particles ma; private double t = 0.0; private static double dt = 5.0*1.0e-15; private static double xMax = 20.0*1.0e-9; private static double yMax = 20.0*1.0e-9; private static double zMax = 20.0*1.0e-9; public static int Nmt = 20*20*20; private int Nsx = 20; private int Nsy = 20; private int Nsz = 20; private int MAT1 = 21; // all private int MAT2 = 20; // Hg private int mat2content = 50; // 0-100 private double rc = ma.getrCutoff(); private double rc2 = rc*rc; private int tempMode = 0; private double contTemp = 1000.0; private int NN = 10000; private int kind[] = new int[NN]; private double xx[] = new double[NN]; private double yy[] = new double[NN]; private double zz[] = new double[NN]; private double vx[] = new double[NN]; private double vy[] = new double[NN]; private double vz[] = new double[NN]; private double ffx[] = new double[NN]; private double ffy[] = new double[NN]; private double ffz[] = new double[NN]; private int reg[][] = new int[NN][300]; private int section[][][][] = new int[Nsx][Nsy][Nsz][100]; private double virial[][] = new double[3][256]; private double memTemp[] = new double[256]; private int virialp = 0; private int resetSW = 0; private int startSW = 1; private int stepSW = 0; private int sleepTime = 1; private int stopSleepTime = 100; private int count = 0; // ------------------------------- class constructor ----------- MolecularDynamics3D(Particles ma2d) { ma = ma2d; setInitialCondition(); } // -------------------------------------- thread run ----------- public void run() { while (true) { count += 1; timeEvolution(); if (sleepTime>0) { try { Thread.sleep(sleepTime); } catch (InterruptedException e) { } } } } // --------------------------- set initial condition ----------- private void setInitialCondition() { int i; double s,mr; count = 0; t = 0.0; s = xMax/20.0; for (i=0; i=ma.getKindMax()) { kind[i] = i%21; } else { kind[i] = MAT1; } xx[i] = 0.8*s*(double)((i%400)%20)+0.5*s; yy[i] = s*(double)((i%400)/20)+0.5*s; zz[i] = s*(double)(i/400)+0.5*s; mr = Math.sqrt(ma.getMassOf(kind[0])/ma.getMassOf(kind[i])); vx[i] = 800.0*mr*(Math.random()-0.5); vy[i] = 800.0*mr*(Math.random()-0.5); vz[i] = 800.0*mr*(Math.random()-0.5); ffx[i] = 0.0; ffy[i] = 0.0; ffz[i] = 0.0; } vAjustment(); initVirial(); } // ---------------------------------- time evolution ----------- private void timeEvolution(){ if (resetSW==1) { setInitialCondition(); resetSW = 0; } if (startSW==1) { timeStep(); } else { if (stepSW==1) { timeStep(); stepSW = 0; } try { Thread.sleep(stopSleepTime); } catch (InterruptedException e) { } } } private void timeStep() { int i; if (tempMode==1) { vAjustment(); } registration(); for (i=0; i<20; i++) { time1Step(); } setVirial(); } private void time1Step() { int i; double dtv2, a2; dtv2 = dt/2.0; t = t + dt; for (i=0; i0.5*xMax) { xj += xMax; }; if (xi-xj<-0.5*xMax) { xj -= xMax; }; if (yi-yj>0.5*yMax) { yj += yMax; }; if (yi-yj<-0.5*yMax) { yj -= yMax; }; if (zi-zj>0.5*zMax) { zj += zMax; }; if (zi-zj<-0.5*zMax) { zj -= zMax; }; r2 = (xi-xj)*(xi-xj)+(yi-yj)*(yi-yj)+(zi-zj)*(zi-zj); if (r21.2){ rr=1.2; }; } for (i=0; i xMax) { xx[i] -= xMax; } if (yy[i] < 0.0) { yy[i] += yMax; } if (yy[i] > yMax) { yy[i] -= yMax; } if (zz[i] < 0.0) { zz[i] += zMax; } if (zz[i] > zMax) { zz[i] -= zMax; } } } //registration private void registration() { int i,j,k,ii,jj,kk,ip,jp,kp,i0,i1,j0,j1,k0,k1,iq; double r2,rreg,rreg2; double xi,xj,yi,yj,zi,zj; preRegistration(); rreg = rc + 3000.0*20*dt; rreg2 = rreg*rreg; for(ip=0;ipip) { xi = xx[ip]; xj = xx[jp]; yi = yy[ip]; yj = yy[jp]; zi = zz[ip]; zj = zz[jp]; if (xi-xj>0.5*xMax) xj += xMax; if (xi-xj<-0.5*xMax) xj -= xMax; if (yi-yj>0.5*yMax) yj += yMax; if (yi-yj<-0.5*yMax) yj -= yMax; if (zi-zj>0.5*zMax) zj += zMax; if (zi-zj<-0.5*zMax) zj -= zMax; r2=(xi-xj)*(xi-xj)+(yi-yj)*(yi-yj)+(zi-zj)*(zi-zj); if (r2=Nsx) i = Nsx-1; j = (int)(Nsy*yy[ip]/yMax); if (j>=Nsy) j = Nsy-1; k = (int)(Nsz*zz[ip]/zMax); if (k>=Nsz) k = Nsz-1; iq = section[i][j][k][0]+1; section[i][j][k][0] = iq; section[i][j][k][iq] = ip; } } // ----------------------------------------- utility ----------- private void vAjustment() { int i; double tmp, r; tmp = temperature(); r = Math.sqrt(contTemp/tmp); for (i=0; im) { m = reg[i][0]; } } return(m); } public int secMax() { int i,j,k,m; m = 0; for(i=0;im) { m = section[i][j][k][0]; } } } } return(m); } // pressure private void initVirial() { int i; for (i=0; i<256; i++) { virial[0][i] = 0.0; virial[1][i] = 0.0; virial[2][i] = 0.0; memTemp[i] = 0.0; } virialp = 0; } private void setVirial() { int i; double xVir,yVir,zVir; xVir = 0.0; yVir = 0.0; zVir = 0.0; for (i=0; i