add the JGF MolDyn benchmark
authorjzhou <jzhou>
Tue, 11 Nov 2008 23:02:32 +0000 (23:02 +0000)
committerjzhou <jzhou>
Tue, 11 Nov 2008 23:02:32 +0000 (23:02 +0000)
Robust/src/Benchmarks/Scheduling/JGFMolDyn/JGFMolDynBench.java [new file with mode: 0644]
Robust/src/Benchmarks/Scheduling/JGFMolDyn/MD.java [new file with mode: 0644]
Robust/src/Benchmarks/Scheduling/JGFMolDyn/MDRunner.java [new file with mode: 0644]
Robust/src/Benchmarks/Scheduling/JGFMolDyn/MyRandom.java [new file with mode: 0644]
Robust/src/Benchmarks/Scheduling/JGFMolDyn/Particle.java [new file with mode: 0644]

diff --git a/Robust/src/Benchmarks/Scheduling/JGFMolDyn/JGFMolDynBench.java b/Robust/src/Benchmarks/Scheduling/JGFMolDyn/JGFMolDynBench.java
new file mode 100644 (file)
index 0000000..c8287f0
--- /dev/null
@@ -0,0 +1,121 @@
+/** Banboo Version */
+
+/**************************************************************************
+*                                                                         *
+*         Java Grande Forum Benchmark Suite - Thread Version 1.0          *
+*                                                                         *
+*                            produced by                                  *
+*                                                                         *
+*                  Java Grande Benchmarking Project                       *
+*                                                                         *
+*                                at                                       *
+*                                                                         *
+*                Edinburgh Parallel Computing Centre                      *
+*                                                                         * 
+*                email: epcc-javagrande@epcc.ed.ac.uk                     *
+*                                                                         *
+*                                                                         *
+*      This version copyright (c) The University of Edinburgh, 2001.      *
+*                         All rights reserved.                            *
+*                                                                         *
+**************************************************************************/
+
+task t1(StartupObject s{initialstate}) {
+    //System.printString("task t1\n");
+    
+    int datasize = 2;
+    int group = 16;
+    
+    MD md = new MD(datasize, group){initialise};
+    
+    taskexit(s{!initialstate});
+}
+
+task t2(MD md{initialise}) {
+    //System.printString("task t2\n");
+    md.initialise();
+    taskexit(md{!initialise, move});
+}
+
+task t3(MD md{move}) {
+    //System.printString("task t3\n");
+    
+    md.domove();
+    md.init();
+    //System.printI(0xd0);
+    for(int i = 0; i < md.group; ++i) {
+       MDRunner runner = new MDRunner(i, md){run};
+       runner.init();
+       //System.printI(0xd1);
+    }
+    //System.printI(0xd2);
+    taskexit(md{!move,update});
+}
+
+/*task t4(MD md{fire}, MDRunner runner{wait}) {
+    System.printString("task t4\n");
+    
+    runner.init();
+    md.counter++;
+    //System.printString("counter: " + md.counter + "\n");
+    if(md.counter == md.group) {
+       //System.printString("Fire finished: " + md.move + "\n");
+       taskexit(md{!fire, update}, runner{!wait, run});
+    } else {
+       taskexit(md{fire}, runner{!wait, run});
+    }
+}*/
+
+task t5(MDRunner runner{run}) {
+    //System.printString("task t5\n");
+    
+    runner.run();
+    
+    taskexit(runner{update, !run});
+}
+
+task t6(MD md{update}, MDRunner runner{update}) {
+    //System.printString("task t6\n");
+    
+    md.update(runner);
+    md.counter++;
+    if(md.counter == md.group) {
+       md.sum();
+       md.scale();
+       if(md.finish()) {
+           taskexit(md{!update, validate}, runner{!update/*, scale*/});
+       } else {
+           taskexit(md{!update, move}, runner{!update/*, scale*/});
+       }
+    } else {
+       taskexit(md{update}, runner{!update/*, scale*/});
+    }
+}
+
+/*task t7(MD md{scale}, MDRunner runner{scale}) {
+    //System.printString("task t7\n");
+    
+    md.counter--;
+    
+    if(md.counter == 0) {
+       if(md.finish()) {
+           // finished
+           taskexit(md{!scale, validate}, runner{!scale, !wait});
+       }
+       taskexit(md{!scale, move}, runner{!scale, !wait});
+    } else {
+       if(md.finish()) {
+           taskexit(md{scale}, runner{!scale, !wait});
+       }
+       taskexit(md{scale}, runner{!scale, !wait});
+    }
+}*/
+
+task t8(MD md{validate}) {
+    //System.printString("task t8\n");
+    
+    md.validate();
+    
+    taskexit(md{!validate, finish});
+}
diff --git a/Robust/src/Benchmarks/Scheduling/JGFMolDyn/MD.java b/Robust/src/Benchmarks/Scheduling/JGFMolDyn/MD.java
new file mode 100644 (file)
index 0000000..75fd7de
--- /dev/null
@@ -0,0 +1,356 @@
+/** Banboo Version */
+
+/**************************************************************************
+ *                                                                         *
+ *             Java Grande Forum Benchmark Suite - Version 2.0             *
+ *                                                                         *
+ *                            produced by                                  *
+ *                                                                         *
+ *                  Java Grande Benchmarking Project                       *
+ *                                                                         *
+ *                                at                                       *
+ *                                                                         *
+ *                Edinburgh Parallel Computing Centre                      *
+ *                                                                         * 
+ *                email: epcc-javagrande@epcc.ed.ac.uk                     *
+ *                                                                         *
+ *                  Original version of this code by                       *
+ *                         Dieter Heermann                                 * 
+ *                       converted to Java by                              *
+ *                Lorna Smith  (l.smith@epcc.ed.ac.uk)                     *
+ *                   (see copyright notice below)                          *
+ *                                                                         *
+ *      This version copyright (c) The University of Edinburgh, 2001.      *
+ *                         All rights reserved.                            *
+ *                                                                         *
+ **************************************************************************/
+
+public class MD {
+    flag initialise;
+    flag move;
+    flag fire;
+    flag update;
+    flag scale;
+    flag validate;
+    flag finish;
+    
+    public int mm;
+    public int group;
+    public int mdsize;
+    public int move;
+    public int movemx;
+    public float side;
+    float hsq, hsq2, den, h, vaver,vaverh;
+    float tref, tscale;
+    int irep;
+    int istop;
+    int iprint;
+    
+    public Particle[] one;
+    public float[][] sh_force;
+    //public float[] epot;
+    //public float[] vir;
+    //public float[] ek;
+    public float epot;
+    public float vir;
+    public float ek;
+    float ekin;
+    
+    public int counter;
+    
+    public void MD(int d, int g) {
+       this.mm = d;
+       this.group = g;
+       this.mdsize = this.mm * this.mm * this.mm * 4;
+       this.one = new Particle [this.mdsize];
+       this.sh_force = new float[3][this.mdsize];
+       this.move = 0;
+       this.movemx = 2;
+       this.den = (float)0.83134;
+       this.h = (float)0.064;
+       this.side = Math.powf((float)(this.mdsize/this.den),(float)0.3333333);
+       this.hsq = this.h * this.h;
+       this.hsq2 = this.hsq * (float)0.5;
+       this.vaver = (float)1.13 * Math.sqrtf((float)this.tref / (float)24.0);
+       this.vaverh = this.vaver * this.h;
+       this.irep = 10;
+       this.istop = 19;
+       //this.epot = (float)0.0;
+       //this.vir = (float)0.0;
+       this.tref = (float)0.722;
+       this.tscale = (float)16.0 / ((float)1.0 * (float)this.mdsize - (float)1.0);
+       this.iprint = 10;
+       this.counter = 0;
+       
+       this.epot = (float)0.0;//new float[this.group + 1];
+       this.vir = (float)0.0;//new float[this.group + 1];
+       this.ek = (float)0.0;//new float[this.group + 1];
+    }
+    
+    public void init() {
+       /*for(int i = 0; i < this.group + 1; i++) {
+           this.epot[i] = (float)0.0;
+           this.vir[i] = (float)0.0;
+       }*/
+       this.epot = (float)0.0;
+       this.vir = (float)0.0;
+       for(int j = 0; j < 3; j++) {
+           for (int i = 0; i < this.mdsize; i++) {
+               this.sh_force[j][i] = (float)0.0;
+           }
+       }
+       this.counter = 0;
+    }
+
+    public void initialise() {
+       
+       /* Particle Generation */
+       float xvelocity, yvelocity, zvelocity;
+       int ijk, lg, i, j, k;
+       float a = this.side/this.mm;
+       xvelocity = (float)0.0;
+       yvelocity = (float)0.0;
+       zvelocity = (float)0.0;
+       //System.printString("here 1\n");
+       //System.printI(0xa0);
+       ijk = 0;
+       for (lg=0; lg<=1; lg++) {
+           for (i=0; i<mm; i++) {
+               for (j=0; j<mm; j++) {
+                   for (k=0; k<mm; k++) {
+                       one[ijk] = new Particle((float)(i*a+lg*a*(float)0.5),(float)(j*a+lg*a*(float)0.5),(k*a),
+                               xvelocity,yvelocity,zvelocity,this);
+                       ijk = ijk + 1;
+                   }
+               }
+           }
+       }
+       for (lg=1; lg<=2; lg++) {
+           for (i=0; i<mm; i++) {
+               for (j=0; j<mm; j++) {
+                   for (k=0; k<mm; k++) {
+                       one[ijk] = new Particle((float)(i*a+(2-lg)*a*(float)0.5),(float)(j*a+(lg-1)*a*(float)0.5),
+                               (float)(k*a+a*(float)0.5),xvelocity,yvelocity,zvelocity,this);
+                       ijk = ijk + 1;
+                   }
+               }
+           }
+       }
+       //System.printString("here 2\n");
+       //System.printI(0xa1);
+       /* Initialise velocities */
+       int iseed;
+       float v1,v2, r;
+       iseed = 0;
+       v1 = (float)0.0;
+       v2 = (float)0.0;
+       
+       MyRandom random = new MyRandom(iseed,v1,v2);
+       //System.printString("here 3\n");
+       //System.printI(0xa2);
+       for (i=0; i<this.mdsize; i+=2) {
+           r  = random.seed();
+           one[i].xvelocity = r*random.v1;
+           one[i+1].xvelocity  = r*random.v2;
+       }
+       //System.printString("here 4\n");
+       //System.printI(0xa3);
+       for (i=0; i<this.mdsize; i+=2) {
+           r  = random.seed();
+           one[i].yvelocity = r*random.v1;
+           one[i+1].yvelocity  = r*random.v2;
+       }
+       //System.printString("here 5\n");
+       //System.printI(0xa4);
+       for (i=0; i<this.mdsize; i+=2) {
+           r  = random.seed();
+           one[i].zvelocity = r*random.v1;
+           one[i+1].zvelocity  = r*random.v2;
+       }
+
+       /*for(i = 0; i < this.mdsize; i++) {
+               System.printString("xvel: " + (int)(one[i].xvelocity*100000) + "; yvel: " + (int)(one[i].yvelocity*100000) + "; zvel: " + (int)(one[i].zvelocity*100000) + "\n");
+       }*/
+       
+       //System.printString("here 6\n");
+       //System.printI(0xa5);
+       /* velocity scaling */
+       float sp, ts, sc;
+       ekin = (float)0.0;
+       sp = (float)0.0;
+       //System.printString("here 7\n");
+       //System.printI(0xa6);
+       for(i=0;i<this.mdsize;i++) {
+           sp = sp + one[i].xvelocity;
+       }
+       sp = sp / this.mdsize;
+       //System.printString("here 8\n");
+       //System.printI(0xa7);
+       for(i=0;i<this.mdsize;i++) {
+           one[i].xvelocity = one[i].xvelocity - sp;
+           ekin = ekin + one[i].xvelocity*one[i].xvelocity;
+       }
+       //System.printString("here 9\n");
+       //System.printI(0xa8);
+       sp = (float)0.0;
+       for(i=0;i<this.mdsize;i++) {
+           sp = sp + one[i].yvelocity;
+       }
+       sp = sp / this.mdsize;
+       //System.printString("here 10\n");
+       //System.printI(0xa9);
+       for(i=0;i<this.mdsize;i++) {
+           one[i].yvelocity = one[i].yvelocity - sp;
+           ekin = ekin + one[i].yvelocity*one[i].yvelocity;
+       }
+       //System.printString("here 11\n");
+       //System.printI(0xa10);
+       sp = (float)0.0;
+       for(i=0;i<this.mdsize;i++) {
+           sp = sp + one[i].zvelocity;
+       }
+       sp = sp / this.mdsize;
+       //System.printString("here 12\n");
+       //System.printI(0xa11);
+       for(i=0;i<this.mdsize;i++) {
+           one[i].zvelocity = one[i].zvelocity - sp;
+           ekin = ekin + one[i].zvelocity*one[i].zvelocity;
+       }
+       //System.printString("here 13\n");
+       //System.printI(0xa12);
+       ts = this.tscale * ekin;
+       sc = h * Math.sqrtf(this.tref/ts);
+
+       for(i=0;i<this.mdsize;i++) {
+           one[i].xvelocity = one[i].xvelocity * sc;     
+           one[i].yvelocity = one[i].yvelocity * sc;     
+           one[i].zvelocity = one[i].zvelocity * sc;     
+       }
+       //System.printString("here 14\n");
+       //System.printI(0xa13);
+    }
+
+
+    public void domove(){
+       for (int i=0;i<this.mdsize;i++) {
+           one[i].domove(this.side,i);       
+       }
+       
+       /*for(int j=0;j<3;j++) {
+           for (int i=0;i<mdsize;i++) {
+               sh_force[j][i] = (float)0.0;
+           }
+       }*/
+       
+       this.move++;
+    }
+    
+    public boolean finish() {
+       if(this.move == this.movemx) {
+           return true;
+       } 
+       return false;
+    }
+
+    public void update(MDRunner runner) {
+       float sum, vel,velt, count, sc;
+       float etot,temp,pres,rp;
+       
+       /* update force arrays */
+       for(int k=0;k<3;k++) {
+           for(int i=0;i<this.mdsize;i++) {
+               sh_force[k][i] += runner.sh_force2[k][i];
+           }
+       }
+       
+       //runner.init();
+
+       //this.epot[runner.id + 1] = runner.epot;
+       //this.vir[runner.id + 1] = runner.vir;
+       this.epot += runner.epot;
+       this.vir += runner.vir;
+    }
+    
+    public void sum() {
+       float sum = (float)0.0; 
+       
+       for (int j=0;j<3;j++) {
+           for (int i=0;i<this.mdsize;i++) {
+               sh_force[j][i] = sh_force[j][i] * hsq2;
+           }
+       }
+       
+       sum = (float)0.0;
+       for (int i=0;i<this.mdsize;i++) {
+           sum = sum + this.one[i].mkekin(hsq2,i);  
+       }
+       ekin = (float)(sum/hsq);
+    }
+    
+    public void scale() {
+       float sum, vel,velt, count, sc;
+       float etot,temp,pres,rp;
+       
+       //runner.epot = this.epot[0];
+       //runner.vir = this.vir[0];
+
+       vel = (float)0.0;
+       count = (float)0.0;
+
+       /* average velocity */
+       for (int i=0;i<this.mdsize;i++) {
+           velt = this.one[i].velavg(vaverh,h);
+           if(velt > vaverh) { 
+               count = count + (float)1.0; 
+           }
+           vel = vel + velt;                    
+       }
+       vel = (float)(vel / h);
+
+       /* temperature scale if required */
+       if((this.move < this.istop) && (((this.move+1) % this.irep) == 0)) {
+           sc = Math.sqrtf(this.tref / (this.tscale*ekin));
+           for (int i=0;i<mdsize;i++) {
+               one[i].dscal(sc,1);
+           }
+           ekin = (float)(this.tref / this.tscale);
+       }
+
+       /* sum to get full potential energy and virial */
+       if(((this.move+1) % this.iprint) == 0) {
+           /*this.ek[runner.id+1] = (float)24.0*ekin;
+           this.epot[runner.id+1] = (float)4.0*this.epot[runner.id+1];
+           etot = this.ek[runner.id+1] + this.epot[runner.id+1];
+           temp = this.tscale * ekin;
+           pres = den * (float)16.0 * (ekin - this.vir[runner.id+1]) / mdsize;
+           vel = vel / this.mdsize; 
+           rp = (count / this.mdsize) * (float)100.0;
+           
+           if(this.counter == this.group) {*/
+               this.ek = (float)24.0*ekin;
+               //this.epot[0] = (float)4.0*this.epot[0];
+               //etot = this.ek[0] + this.epot[0];
+               //temp = this.tscale * ekin;
+               //pres = den * (float)16.0 * (ekin - this.vir[0]) / mdsize;
+               //vel = vel / this.mdsize; 
+               //rp = (count / this.mdsize) * (float)100.0;
+               //System.printString("ek: " + (int)(this.ek*1000000) + " (" + this.move + ")\n");
+           //}
+       }
+    }
+    
+    public void validate() {
+       float refval = (float)1731.4306625334357;
+       float dev = Math.abs(this.ek - refval);
+       if (dev > 1.0e-10 ){
+           //System.printString("Validation failed\n");
+           //System.printString("Kinetic Energy = " + (int)(this.ek*1000000) + ";  " + (int)(refval*1000000) + "\n");
+           //System.printI(0xdddf);
+           //System.printI((int)(this.ek[0]*1000000));
+           //System.printI((int)(refval*1000000));
+           //System.printI((int)(ek*10000));
+           //System.printI((int)(refval*10000));
+       }
+    }
+}
diff --git a/Robust/src/Benchmarks/Scheduling/JGFMolDyn/MDRunner.java b/Robust/src/Benchmarks/Scheduling/JGFMolDyn/MDRunner.java
new file mode 100644 (file)
index 0000000..7b75551
--- /dev/null
@@ -0,0 +1,72 @@
+public class MDRunner {
+    flag wait;
+    flag run;
+    flag update;
+    flag scale;
+    
+    public int id;  
+    int mm;
+    int mdsize;
+    public int group;
+    public int ilow, iupper;
+    //float l;
+    float rcoff;
+    //rcoffs,
+    float side/*sideh,*/;
+    //float rand;
+    public float [][] sh_force2;
+    public float epot;
+    public float vir;
+    public float ek;
+    //int npartm,tint;
+    //int iprint = 10;
+    MD md;
+    Particle[] one;
+
+    public MDRunner(int id, MD m) {
+       this.id=id;
+       this.md = m;
+       this.mm=this.md.mm;
+       //this.l = 50e-10;
+       this.mdsize = this.md.mdsize;
+       this.group = this.md.group;
+       this.side = this.md.side;
+       this.rcoff = this.mm/(float)4.0;
+//     this.rcoffs = this.rcoff * this.rcoff;
+       //this.sideh = this.side * 0.5;
+       //this.npartm = this.mdsize - 1;
+       //this.iprint = 10;
+       int slice = (this.mdsize - 1) / this.group + 1; 
+       this.ilow = this.id * slice; 
+       this.iupper = (this.id+1) * slice; 
+       if (this.iupper > this.mdsize )  {
+           iupper = this.mdsize; 
+       }
+       sh_force2 = new float[3][this.mdsize];
+       
+       this.one = this.md.one;
+       this.epot = (float)0.0;//this.md.epot[id+1];
+       this.vir = (float)0.0;//this.md.vir[id+1];
+       this.ek = (float)0.0;//this.md.ek[id+1];
+    }
+    
+    public void init() {
+       this.epot = (float)0.0;
+       this.vir = (float)0.0;
+       for(int i = 0; i < this.mdsize; i++) {
+           for(int j = 0; j < 3; j++) {
+               this.sh_force2[j][i] = (float)0.0;
+           }
+       }
+    }
+
+    public void run() {
+       /* compute forces */
+       //System.printString("here 1: " + this.id + "\n");
+       for (int i = this.ilow; i < this.iupper; i++) {
+           one[i].force(side,rcoff,mdsize,i,this); 
+       }
+       //System.printString("here 2: " + this.id + "\n");
+    }
+
+}
\ No newline at end of file
diff --git a/Robust/src/Benchmarks/Scheduling/JGFMolDyn/MyRandom.java b/Robust/src/Benchmarks/Scheduling/JGFMolDyn/MyRandom.java
new file mode 100644 (file)
index 0000000..3fdb7e4
--- /dev/null
@@ -0,0 +1,61 @@
+public class MyRandom {
+
+    public int iseed;
+    public float v1,v2;
+
+    public MyRandom(int iseed, float v1, float v2) {
+       this.iseed = iseed;
+       this.v1 = v1;
+       this.v2 = v2;
+    }
+
+    public float update() {
+
+       float rand;
+       float scale= (float)4.656612875e-10;
+
+       int is1,is2,iss2;
+       int imult= 16807;
+       int imod = 2147483647;
+
+       if (iseed<=0) { 
+           iseed = 1; 
+           }
+
+       is2 = iseed % 32768;
+       is1 = (iseed-is2)/32768;
+       iss2 = is2 * imult;
+       is2 = iss2 % 32768;
+       is1 = (is1 * imult+(iss2-is2) / 32768) % (65536);
+
+       iseed = (is1 * 32768 + is2) % imod;
+
+       rand = scale * iseed;
+
+       return rand;
+
+    }
+
+    public float seed() {
+
+       float s,u1,u2,r;
+       s = (float)1.0;
+       //do {
+           //System.printI(0xb1);
+           u1 = update();
+           u2 = update();
+
+       //System.printI(0xb2);
+           v1 = (float)2.0 * u1 - (float)1.0;
+           v2 = (float)2.0 * u2 - (float)1.0;
+           s = v1*v1 + v2*v2;
+       //System.printI(0xb3);
+       //} while (s >= (float)1.0);
+        s = s - (int)s;
+       //System.printI(0xb4);
+       r = Math.sqrtf((float)(-2.0*Math.logf(s))/(float)s);
+       //System.printI(0xb5);
+       return r;
+
+    }
+}
diff --git a/Robust/src/Benchmarks/Scheduling/JGFMolDyn/Particle.java b/Robust/src/Benchmarks/Scheduling/JGFMolDyn/Particle.java
new file mode 100644 (file)
index 0000000..5376afa
--- /dev/null
@@ -0,0 +1,139 @@
+public class Particle {
+
+    public float xcoord, ycoord, zcoord;
+    public float xvelocity,yvelocity,zvelocity;
+    int id;
+    //float [][] sh_force;
+    //float [][][] sh_force2;
+    MD md;
+
+    public Particle(float xcoord, float ycoord, float zcoord, float xvelocity,
+           float yvelocity,float zvelocity,/*float [][] sh_force, 
+           float [][][] sh_force2, int id, */MD m) {
+
+       this.xcoord = xcoord; 
+       this.ycoord = ycoord; 
+       this.zcoord = zcoord;
+       this.xvelocity = xvelocity;
+       this.yvelocity = yvelocity;
+       this.zvelocity = zvelocity;
+       //this.sh_force = sh_force;
+       //this.sh_force2 = sh_force2;
+       this.id=id;
+       this.md=m;
+    }
+
+    public void domove(float side,int part_id) {
+       xcoord = xcoord + xvelocity + this.md.sh_force[0][part_id];
+       ycoord = ycoord + yvelocity + this.md.sh_force[1][part_id];
+       zcoord = zcoord + zvelocity + this.md.sh_force[2][part_id];
+
+       if(xcoord < 0) { xcoord = xcoord + side; } 
+       if(xcoord > side) { xcoord = xcoord - side; }
+       if(ycoord < 0) { ycoord = ycoord + side; }
+       if(ycoord > side) { ycoord = ycoord - side; }
+       if(zcoord < 0) { zcoord = zcoord + side; }
+       if(zcoord > side) { zcoord = zcoord - side; }
+
+       xvelocity = xvelocity + this.md.sh_force[0][part_id];
+       yvelocity = yvelocity + this.md.sh_force[1][part_id];
+       zvelocity = zvelocity + this.md.sh_force[2][part_id];
+       //System.printI(0xc0);
+    }
+
+    public void force(float side, float rcoff,int mdsize,int x, MDRunner runner) {
+       float sideh;
+       float rcoffs;
+
+       float fxi,fyi,fzi;
+       float rd,rrd,rrd2,rrd3,rrd4,rrd6,rrd7,r148;
+       float forcex,forcey,forcez;
+       
+       float xx, yy, zz;
+
+       sideh = (float)0.5*side; 
+       rcoffs = rcoff*rcoff;
+
+       fxi = (float)0.0;
+       fyi = (float)0.0;
+       fzi = (float)0.0;
+       //System.printString("here 111: " + runner.id + "\n");
+       for (int i=x+1;i<mdsize;i++) {
+           xx = this.xcoord - this.md.one[i].xcoord;
+           yy = this.ycoord - this.md.one[i].ycoord;
+           zz = this.zcoord - this.md.one[i].zcoord;
+
+           if(xx < (-sideh)) { xx = xx + side; }
+           if(xx > (sideh))  { xx = xx - side; }
+           if(yy < (-sideh)) { yy = yy + side; }
+           if(yy > (sideh))  { yy = yy - side; }
+           if(zz < (-sideh)) { zz = zz + side; }
+           if(zz > (sideh))  { zz = zz - side; }
+
+
+           rd = xx*xx + yy*yy + zz*zz;
+
+           if(rd <= rcoffs) {
+               rrd = (float)1.0/rd;
+               rrd2 = rrd*rrd;
+               rrd3 = rrd2*rrd;
+               rrd4 = rrd2*rrd2;
+               rrd6 = rrd2*rrd4;
+               rrd7 = rrd6*rrd;
+               runner.epot = runner.epot + (rrd6 - rrd3);
+               r148 = rrd7 - (float)0.5*rrd4;
+               runner.vir = runner.vir - rd*r148;
+               forcex = xx * r148;
+               fxi = fxi + forcex;
+
+               runner.sh_force2[0][i] = runner.sh_force2[0][i] - forcex;
+
+               forcey = yy * r148;
+               fyi = fyi + forcey;
+
+               runner.sh_force2[1][i] = runner.sh_force2[1][i] - forcey;
+
+               forcez = zz * r148;
+               fzi = fzi + forcez;
+
+               runner.sh_force2[2][i] = runner.sh_force2[2][i] - forcez;
+
+               //this.md.interacts[id]++;
+           }
+
+       }
+       //System.printString("here 222: " + runner.id + "\n");
+       runner.sh_force2[0][x] = runner.sh_force2[0][x] + fxi;
+       runner.sh_force2[1][x] = runner.sh_force2[1][x] + fyi;
+       runner.sh_force2[2][x] = runner.sh_force2[2][x] + fzi;
+       //System.printString("here 333: " + runner.id + "\n");
+    }
+
+    public float mkekin(float hsq2,int part_id) {
+       float sumt = (float)0.0; 
+
+       xvelocity = xvelocity + this.md.sh_force[0][part_id]; 
+       yvelocity = yvelocity + this.md.sh_force[1][part_id]; 
+       zvelocity = zvelocity + this.md.sh_force[2][part_id]; 
+
+       sumt = (xvelocity*xvelocity)+(yvelocity*yvelocity)+(zvelocity*zvelocity);
+       return sumt;
+    }
+
+    public float velavg(float vaverh,float h) {
+       float velt;
+       float sq;
+
+       sq = Math.sqrtf(xvelocity*xvelocity + yvelocity*yvelocity + zvelocity*zvelocity);
+
+       velt = sq;
+       return velt;
+    }
+
+    public void dscal(float sc,int incx) {
+       xvelocity = xvelocity * sc;
+       yvelocity = yvelocity * sc;   
+       zvelocity = zvelocity * sc;   
+    }
+
+}
\ No newline at end of file