ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/SimInfo.cpp
(Generate patch)

Comparing trunk/OOPSE/libmdtools/SimInfo.cpp (file contents):
Revision 588 by gezelter, Thu Jul 10 17:10:56 2003 UTC vs.
Revision 621 by gezelter, Wed Jul 16 02:11:02 2003 UTC

# Line 34 | Line 34 | SimInfo::SimInfo(){
34    setTemp = 0;
35    thermalTime = 0.0;
36    rCut = 0.0;
37 +  ecr = 0.0;
38 +  est = 0.0;
39  
40    usePBC = 0;
41    useLJ = 0;
# Line 77 | Line 79 | void SimInfo::setBoxM( double theBox[3][3] ){
79    for(i=0; i < 3; i++)
80      for (j=0; j < 3; j++) Hmat[i][j] = theBox[i][j];
81    
82 <  cerr
83 <    << "setting Hmat ->\n"
84 <    << "[ " << Hmat[0][0] << ", " << Hmat[0][1] << ", " << Hmat[0][2] << " ]\n"
85 <    << "[ " << Hmat[1][0] << ", " << Hmat[1][1] << ", " << Hmat[1][2] << " ]\n"
86 <    << "[ " << Hmat[2][0] << ", " << Hmat[2][1] << ", " << Hmat[2][2] << " ]\n";
82 >  //  cerr
83 >  // << "setting Hmat ->\n"
84 >  // << "[ " << Hmat[0][0] << ", " << Hmat[0][1] << ", " << Hmat[0][2] << " ]\n"
85 >  // << "[ " << Hmat[1][0] << ", " << Hmat[1][1] << ", " << Hmat[1][2] << " ]\n"
86 >  // << "[ " << Hmat[2][0] << ", " << Hmat[2][1] << ", " << Hmat[2][2] << " ]\n";
87  
88    calcBoxL();
89    calcHmatInv();
# Line 93 | Line 95 | void SimInfo::setBoxM( double theBox[3][3] ){
95      }
96    }
97  
98 <  setFortranBoxSize(FortranHmat, FortranHmatI, &orthoRhombic);
98 >  setFortranBoxSize(FortranHmat, FortranHmatInv, &orthoRhombic);
99  
100 <  smallestBoxL = boxLx;
101 <  if (boxLy < smallestBoxL) smallestBoxL = boxLy;
102 <  if (boxLz < smallestBoxL) smallestBoxL = boxLz;
100 >  smallestBoxL = boxL[0];
101 >  if (boxL[1] < smallestBoxL) smallestBoxL = boxL[1];
102 >  if (boxL[2] > smallestBoxL) smallestBoxL = boxL[2];
103  
104    maxCutoff = smallestBoxL / 2.0;
105  
# Line 107 | Line 109 | void SimInfo::setBoxM( double theBox[3][3] ){
109               maxCutoff );
110      painCave.isFatal = 0;
111      simError();
110 –
112      rList = maxCutoff;
113  
114 <    sprintf( painCave.errMsg,
114 <             "New Box size is forcing cutoff radius down to %lf\n",
115 <             maxCutoff - 1.0 );
116 <    painCave.isFatal = 0;
117 <    simError();
118 <
119 <    rCut = rList - 1.0;
120 <
121 <    // list radius changed so we have to refresh the simulation structure.
122 <    refreshSim();
123 <  }
124 <
125 <  if (rCut > maxCutoff) {
126 <    sprintf( painCave.errMsg,
127 <             "New Box size is forcing cutoff radius down to %lf\n",
128 <             maxCutoff );
129 <    painCave.isFatal = 0;
130 <    simError();
131 <
132 <    status = 0;
133 <    LJ_new_rcut(&rCut, &status);
134 <    if (status != 0) {
114 >    if (rCut > (rList - 1.0)) {
115        sprintf( painCave.errMsg,
116 <               "Error in recomputing LJ shifts based on new rcut\n");
117 <      painCave.isFatal = 1;
116 >               "New Box size is forcing LJ cutoff radius down to %lf\n",
117 >               rList - 1.0 );
118 >      painCave.isFatal = 0;
119        simError();
120 +      rCut = rList - 1.0;
121      }
122 <  }
122 >
123 >    if( ecr > (rList - 1.0) ){
124 >      sprintf( painCave.errMsg,
125 >               "New Box size is forcing electrostaticCutoffRadius "
126 >               "down to %lf\n"
127 >               "electrostaticSkinThickness is now %lf\n",
128 >               rList - 1.0, 0.05*(rList-1.0) );
129 >      painCave.isFatal = 0;
130 >      simError();      
131 >      ecr = maxCutoff;
132 >      est = 0.05 * ecr;
133 >    }
134 >
135 >    // At least one of the radii changed, so we need a refresh:
136 >    refreshSim();
137 >  }    
138   }
139  
140  
# Line 153 | Line 150 | void SimInfo::scaleBox(double scale) {
150    double theBox[3][3];
151    int i, j;
152  
153 <  cerr << "Scaling box by " << scale << "\n";
153 >  // cerr << "Scaling box by " << scale << "\n";
154  
155    for(i=0; i<3; i++)
156      for (j=0; j<3; j++) theBox[i][j] = Hmat[i][j]*scale;
# Line 163 | Line 160 | void SimInfo::calcHmatInv( void ) {
160   }
161  
162   void SimInfo::calcHmatInv( void ) {
163 <
163 >  
164 >  int i,j;
165    double smallDiag;
166    double tol;
167    double sanity[3][3];
# Line 173 | Line 171 | void SimInfo::calcHmatInv( void ) {
171    // Check the inverse to make sure it is sane:
172  
173    matMul3( Hmat, HmatInv, sanity );
176 –
177 –  cerr << "sanity => \n"
178 –       << sanity[0][0] << "\t" << sanity[0][1] << "\t" << sanity [0][2] << "\n"
179 –       << sanity[1][0] << "\t" << sanity[1][1] << "\t" << sanity [1][2] << "\n"
180 –       << sanity[2][0] << "\t" << sanity[2][1] << "\t" << sanity [2][2]
181 –       << "\n";
174      
175    // check to see if Hmat is orthorhombic
176    
# Line 271 | Line 263 | void SimInfo::matVecMul3(double m[3][3], double inVec[
263    outVec[1] = m[1][0]*a0 + m[1][1]*a1 + m[1][2]*a2;
264    outVec[2] = m[2][0]*a0 + m[2][1]*a1 + m[2][2]*a2;
265   }
266 +
267 + void SimInfo::transposeMat3(double in[3][3], double out[3][3]) {
268 +  double temp[3][3];
269 +  int i, j;
270 +
271 +  for (i = 0; i < 3; i++) {
272 +    for (j = 0; j < 3; j++) {
273 +      temp[j][i] = in[i][j];
274 +    }
275 +  }
276 +  for (i = 0; i < 3; i++) {
277 +    for (j = 0; j < 3; j++) {
278 +      out[i][j] = temp[i][j];
279 +    }
280 +  }
281 + }
282    
283 + void SimInfo::printMat3(double A[3][3] ){
284 +
285 +  std::cerr
286 +            << "[ " << A[0][0] << ", " << A[0][1] << ", " << A[0][2] << " ]\n"
287 +            << "[ " << A[1][0] << ", " << A[1][1] << ", " << A[1][2] << " ]\n"
288 +            << "[ " << A[2][0] << ", " << A[2][1] << ", " << A[2][2] << " ]\n";
289 + }
290 +
291 + void SimInfo::printMat9(double A[9] ){
292 +
293 +  std::cerr
294 +            << "[ " << A[0] << ", " << A[1] << ", " << A[2] << " ]\n"
295 +            << "[ " << A[3] << ", " << A[4] << ", " << A[5] << " ]\n"
296 +            << "[ " << A[6] << ", " << A[7] << ", " << A[8] << " ]\n";
297 + }
298 +
299   void SimInfo::calcBoxL( void ){
300  
301    double dx, dy, dz, dsq;
# Line 285 | Line 309 | void SimInfo::calcBoxL( void ){
309    
310    dx = Hmat[0][0]; dy = Hmat[1][0]; dz = Hmat[2][0];
311    dsq = dx*dx + dy*dy + dz*dz;
312 <  boxLx = sqrt( dsq );
312 >  boxL[0] = sqrt( dsq );
313  
314    // boxLy
315    
316    dx = Hmat[0][1]; dy = Hmat[1][1]; dz = Hmat[2][1];
317    dsq = dx*dx + dy*dy + dz*dz;
318 <  boxLy = sqrt( dsq );
318 >  boxL[1] = sqrt( dsq );
319  
320    // boxLz
321    
322    dx = Hmat[0][2]; dy = Hmat[1][2]; dz = Hmat[2][2];
323    dsq = dx*dx + dy*dy + dz*dz;
324 <  boxLz = sqrt( dsq );
324 >  boxL[2] = sqrt( dsq );
325    
326   }
327  

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines