Skip to content
Snippets Groups Projects
Commit 696c7d87 authored by Gareth Tribello's avatar Gareth Tribello
Browse files

Fixed small bug in pesmd

parent 803a1f4d
No related branches found
No related tags found
No related merge requests found
......@@ -212,7 +212,7 @@ public:
double tke=0;
for(int i=0; i<nat; ++i) {
for(int j=0; j<3; ++j) {
if( 3*i+j>dim ) break;
if( 3*i+j>dim-1 ) break;
tke += 0.5*velocities[i][j]*velocities[i][j];
}
}
......@@ -241,7 +241,7 @@ public:
double lrand=sqrt((1.-lscale*lscale)*temp);
for(int j=0; j<nat; ++j) {
for(int k=0; k<3; ++k) {
if( 3*j+k>dim ) break;
if( 3*j+k>dim-1 ) break;
therm_eng=therm_eng+0.5*velocities[j][k]*velocities[j][k];
velocities[j][k]=lscale*velocities[j][k]+lrand*random.Gaussian();
therm_eng=therm_eng-0.5*velocities[j][k]*velocities[0][k];
......@@ -251,7 +251,7 @@ public:
// First step of velocity verlet
for(int j=0; j<nat; ++j) {
for(int k=0; k<3; ++k) {
if( 3*j+k>dim ) break;
if( 3*j+k>dim-1 ) break;
velocities[j][k] = velocities[j][k] + 0.5*tstep*forces[1+j][k];
positions[1+j][k] = positions[1+j][k] + tstep*velocities[j][k];
}
......@@ -274,7 +274,7 @@ public:
// Second step of velocity verlet
for(int j=0; j<nat; ++j) {
for(int k=0; k<3; ++k) {
if( 3*j+k>dim ) break;
if( 3*j+k>dim-1 ) break;
velocities[j][k] = velocities[j][k] + 0.5*tstep*forces[1+j][k];
}
}
......@@ -284,7 +284,7 @@ public:
lrand=sqrt((1.-lscale*lscale)*temp);
for(int j=0; j<nat; ++j) {
for(int k=0; k<3; ++k) {
if( 3*j+k>dim ) break;
if( 3*j+k>dim-1 ) break;
therm_eng=therm_eng+0.5*velocities[j][k]*velocities[j][k];
velocities[j][k]=lscale*velocities[j][k]+lrand*random.Gaussian();
therm_eng=therm_eng-0.5*velocities[j][k]*velocities[j][k];
......@@ -294,7 +294,7 @@ public:
tke=0;
for(int i=0; i<nat; ++i) {
for(int j=0; j<3; ++j) {
if( 3*i+j>dim ) break;
if( 3*i+j>dim-1 ) break;
tke += 0.5*velocities[i][j]*velocities[i][j];
}
}
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment