| 9 | 
  | 
  thetaStep = PI / bandwidth; | 
| 10 | 
  | 
  thetaMin = thetaStep / 2.0; | 
| 11 | 
  | 
  phiStep = thetaStep * 2.0; | 
| 12 | 
– | 
         | 
| 13 | 
– | 
  //zero out the rot mats | 
| 14 | 
– | 
  for (i=0; i<3; i++) { | 
| 15 | 
– | 
    for (j=0; j<3; j++) { | 
| 16 | 
– | 
      rotX[i][j] = 0.0; | 
| 17 | 
– | 
      rotZ[i][j] = 0.0; | 
| 18 | 
– | 
      rbMatrix[i][j] = 0.0; | 
| 19 | 
– | 
    } | 
| 20 | 
– | 
  } | 
| 12 | 
  | 
} | 
| 13 | 
  | 
 | 
| 14 | 
  | 
GridBuilder::~GridBuilder() { | 
| 15 | 
  | 
} | 
| 16 | 
  | 
 | 
| 17 | 
< | 
void GridBuilder::launchProbe(int forceField, vector<double> sigmaGrid, vector<double> sGrid, | 
| 18 | 
< | 
                              vector<double> epsGrid){ | 
| 17 | 
> | 
void GridBuilder::launchProbe(int forceField, vector<double> sigmaGrid,  | 
| 18 | 
> | 
                              vector<double> sGrid, vector<double> epsGrid){ | 
| 19 | 
  | 
  ofstream sigmaOut("sigma.grid"); | 
| 20 | 
  | 
  ofstream sOut("s.grid"); | 
| 21 | 
  | 
  ofstream epsOut("eps.grid"); | 
| 34 | 
  | 
  if (startDist < minDist) | 
| 35 | 
  | 
    startDist = minDist; | 
| 36 | 
  | 
 | 
| 37 | 
+ | 
  printf("startDist = %lf\n", startDist); | 
| 38 | 
+ | 
 | 
| 39 | 
  | 
  //set the initial orientation of the body and loop over theta values | 
| 40 | 
< | 
  phiVal = 0.0; | 
| 41 | 
< | 
  thetaVal = thetaMin; | 
| 42 | 
< | 
  rotBody(phiVal, thetaVal); | 
| 43 | 
< | 
  for (k=0; k<bandwidth; k++){   | 
| 44 | 
< | 
        //loop over phi values starting with phi = 0.0 | 
| 45 | 
< | 
    for (j=0; j<bandwidth; j++){ | 
| 40 | 
> | 
 | 
| 41 | 
> | 
  for (k =0; k < bandwidth; k++) { | 
| 42 | 
> | 
    thetaVal = thetaMin + k*thetaStep; | 
| 43 | 
> | 
    for (j=0; j < bandwidth; j++) { | 
| 44 | 
> | 
      phiVal = j*phiStep; | 
| 45 | 
> | 
 | 
| 46 | 
> | 
      printf("setting Euler, phi = %lf\ttheta = %lf\n", phiVal, thetaVal); | 
| 47 | 
> | 
 | 
| 48 | 
> | 
      rbMol->setEuler(0.0, thetaVal, phiVal); | 
| 49 | 
> | 
 | 
| 50 | 
  | 
      releaseProbe(startDist); | 
| 51 | 
  | 
 | 
| 52 | 
+ | 
      printf("found sigDist = %lf\t sDist = %lf \t epsVal = %lf\n",  | 
| 53 | 
+ | 
             sigDist, sDist, epsVal); | 
| 54 | 
+ | 
 | 
| 55 | 
  | 
      sigList.push_back(sigDist); | 
| 56 | 
  | 
      sList.push_back(sDist); | 
| 57 | 
  | 
      epsList.push_back(epsVal); | 
| 58 | 
< | 
       | 
| 59 | 
< | 
      phiVal += phiStep; | 
| 60 | 
< | 
      rotBody(phiVal, thetaVal); | 
| 58 | 
> | 
 | 
| 59 | 
  | 
    } | 
| 62 | 
– | 
    phiVal = 0.0; | 
| 63 | 
– | 
    thetaVal += thetaStep; | 
| 64 | 
– | 
    rotBody(phiVal, thetaVal); | 
| 65 | 
– | 
    printf("step theta %i\n",k); | 
| 60 | 
  | 
  }              | 
| 61 | 
  | 
} | 
| 62 | 
  | 
 | 
| 72 | 
  | 
  tooClose = 0; | 
| 73 | 
  | 
  epsVal = 0; | 
| 74 | 
  | 
  rhoStep = 0.1; //the distance the probe atom moves between steps | 
| 75 | 
< | 
         | 
| 82 | 
< | 
         | 
| 75 | 
> | 
                 | 
| 76 | 
  | 
  while (!tooClose){ | 
| 77 | 
  | 
    calcEnergy(); | 
| 78 | 
  | 
    potProgress.push_back(potEnergy); | 
| 142 | 
  | 
  } | 
| 143 | 
  | 
   | 
| 144 | 
  | 
  potEnergy = 0.0; | 
| 145 | 
+ | 
 | 
| 146 | 
+ | 
  rbMol->getAtomPos(rbAtomPos, 0); | 
| 147 | 
+ | 
 | 
| 148 | 
+ | 
  printf("atom0 pos = %lf\t%lf\t%lf\n", rbAtomPos[0], rbAtomPos[1], rbAtomPos[2]); | 
| 149 | 
+ | 
 | 
| 150 | 
+ | 
 | 
| 151 | 
  | 
   | 
| 152 | 
  | 
  for(i=0; i<rbMol->getNumAtoms(); i++){ | 
| 153 | 
  | 
    rbMol->getAtomPos(rbAtomPos, i); | 
| 215 | 
  | 
  } | 
| 216 | 
  | 
}  | 
| 217 | 
  | 
 | 
| 219 | 
– | 
void GridBuilder::rotBody(double pValue, double tValue){ | 
| 220 | 
– | 
  //zero out the euler angles | 
| 221 | 
– | 
  for (l=0; l<3; l++) | 
| 222 | 
– | 
    angles[i] = 0.0; | 
| 223 | 
– | 
         | 
| 224 | 
– | 
  //the phi euler angle is for rotation about the z-axis (we use the zxz convention) | 
| 225 | 
– | 
  angles[0] = pValue; | 
| 226 | 
– | 
  //the second euler angle is for rotation about the x-axis (we use the zxz convention) | 
| 227 | 
– | 
  angles[1] = tValue; | 
| 228 | 
– | 
         | 
| 229 | 
– | 
  //obtain the rotation matrix through the rigid body class | 
| 230 | 
– | 
  rbMol->doEulerToRotMat(angles, rotX); | 
| 231 | 
– | 
   | 
| 232 | 
– | 
  //start from the reference position | 
| 233 | 
– | 
  identityMat3(rbMatrix); | 
| 234 | 
– | 
  rbMol->setA(rbMatrix); | 
| 235 | 
– | 
   | 
| 236 | 
– | 
  //rotate the rigid body | 
| 237 | 
– | 
  matMul3(rotX, rbMatrix, rotatedMat); | 
| 238 | 
– | 
  rbMol->setA(rotatedMat);       | 
| 239 | 
– | 
} | 
| 240 | 
– | 
 | 
| 218 | 
  | 
void GridBuilder::printGridFiles(){ | 
| 219 | 
  | 
  ofstream sigmaOut("sigma.grid"); | 
| 220 | 
  | 
  ofstream sOut("s.grid"); | 
| 225 | 
  | 
    sOut << sList[k] << "\n0\n";     | 
| 226 | 
  | 
    epsOut << epsList[k] << "\n0\n"; | 
| 227 | 
  | 
  } | 
| 228 | 
< | 
} | 
| 228 | 
> | 
} |