| 1 |
< |
#include <cstdlib> |
| 2 |
< |
#include <cstring> |
| 3 |
< |
#include <cmath> |
| 1 |
> |
#include <stdlib.h> |
| 2 |
> |
#include <string.h> |
| 3 |
> |
#include <math.h> |
| 4 |
|
|
| 5 |
|
#include <iostream> |
| 6 |
|
using namespace std; |
| 94 |
|
|
| 95 |
|
void SimInfo::setBoxM( double theBox[3][3] ){ |
| 96 |
|
|
| 97 |
< |
int i, j, status; |
| 98 |
< |
double smallestBoxL, maxCutoff; |
| 97 |
> |
int i, j; |
| 98 |
|
double FortranHmat[9]; // to preserve compatibility with Fortran the |
| 99 |
|
// ordering in the array is as follows: |
| 100 |
|
// [ 0 3 6 ] |
| 146 |
|
|
| 147 |
|
void SimInfo::calcHmatInv( void ) { |
| 148 |
|
|
| 149 |
+ |
int oldOrtho; |
| 150 |
|
int i,j; |
| 151 |
|
double smallDiag; |
| 152 |
|
double tol; |
| 154 |
|
|
| 155 |
|
invertMat3( Hmat, HmatInv ); |
| 156 |
|
|
| 157 |
– |
// Check the inverse to make sure it is sane: |
| 158 |
– |
|
| 159 |
– |
matMul3( Hmat, HmatInv, sanity ); |
| 160 |
– |
|
| 157 |
|
// check to see if Hmat is orthorhombic |
| 158 |
|
|
| 159 |
< |
smallDiag = Hmat[0][0]; |
| 160 |
< |
if(smallDiag > Hmat[1][1]) smallDiag = Hmat[1][1]; |
| 161 |
< |
if(smallDiag > Hmat[2][2]) smallDiag = Hmat[2][2]; |
| 159 |
> |
oldOrtho = orthoRhombic; |
| 160 |
> |
|
| 161 |
> |
smallDiag = fabs(Hmat[0][0]); |
| 162 |
> |
if(smallDiag > fabs(Hmat[1][1])) smallDiag = fabs(Hmat[1][1]); |
| 163 |
> |
if(smallDiag > fabs(Hmat[2][2])) smallDiag = fabs(Hmat[2][2]); |
| 164 |
|
tol = smallDiag * 1E-6; |
| 165 |
|
|
| 166 |
|
orthoRhombic = 1; |
| 169 |
|
for (j = 0 ; j < 3; j++) { |
| 170 |
|
if (i != j) { |
| 171 |
|
if (orthoRhombic) { |
| 172 |
< |
if (Hmat[i][j] >= tol) orthoRhombic = 0; |
| 172 |
> |
if ( fabs(Hmat[i][j]) >= tol) orthoRhombic = 0; |
| 173 |
|
} |
| 174 |
|
} |
| 175 |
|
} |
| 176 |
|
} |
| 177 |
+ |
|
| 178 |
+ |
if( oldOrtho != orthoRhombic ){ |
| 179 |
+ |
|
| 180 |
+ |
if( orthoRhombic ){ |
| 181 |
+ |
sprintf( painCave.errMsg, |
| 182 |
+ |
"Hmat is switching from Non-Orthorhombic to OrthoRhombic\n" |
| 183 |
+ |
" If this is a bad thing change the ortho tolerance in SimInfo.\n" ); |
| 184 |
+ |
simError(); |
| 185 |
+ |
} |
| 186 |
+ |
else { |
| 187 |
+ |
sprintf( painCave.errMsg, |
| 188 |
+ |
"Hmat is switching from Orthorhombic to Non-OrthoRhombic\n" |
| 189 |
+ |
" If this is a bad thing change the ortho tolerance in SimInfo.\n" ); |
| 190 |
+ |
simError(); |
| 191 |
+ |
} |
| 192 |
+ |
} |
| 193 |
|
} |
| 194 |
|
|
| 195 |
|
double SimInfo::matDet3(double a[3][3]) { |
| 316 |
|
void SimInfo::calcBoxL( void ){ |
| 317 |
|
|
| 318 |
|
double dx, dy, dz, dsq; |
| 305 |
– |
int i; |
| 319 |
|
|
| 320 |
|
// boxVol = Determinant of Hmat |
| 321 |
|
|
| 385 |
|
|
| 386 |
|
void SimInfo::wrapVector( double thePos[3] ){ |
| 387 |
|
|
| 388 |
< |
int i, j, k; |
| 388 |
> |
int i; |
| 389 |
|
double scaled[3]; |
| 390 |
|
|
| 391 |
|
if( !orthoRhombic ){ |
| 423 |
|
|
| 424 |
|
|
| 425 |
|
int SimInfo::getNDF(){ |
| 426 |
< |
int ndf_local, ndf; |
| 426 |
> |
int ndf_local; |
| 427 |
|
|
| 428 |
|
ndf_local = 3 * n_atoms + 3 * n_oriented - n_constraints; |
| 429 |
|
|
| 439 |
|
} |
| 440 |
|
|
| 441 |
|
int SimInfo::getNDFraw() { |
| 442 |
< |
int ndfRaw_local, ndfRaw; |
| 442 |
> |
int ndfRaw_local; |
| 443 |
|
|
| 444 |
|
// Raw degrees of freedom that we have to set |
| 445 |
|
ndfRaw_local = 3 * n_atoms + 3 * n_oriented; |
| 454 |
|
} |
| 455 |
|
|
| 456 |
|
int SimInfo::getNDFtranslational() { |
| 457 |
< |
int ndfTrans_local, ndfTrans; |
| 457 |
> |
int ndfTrans_local; |
| 458 |
|
|
| 459 |
|
ndfTrans_local = 3 * n_atoms - n_constraints; |
| 460 |
|
|
| 530 |
|
|
| 531 |
|
void SimInfo::setRcut( double theRcut ){ |
| 532 |
|
|
| 520 |
– |
if( !haveOrigRcut ){ |
| 521 |
– |
haveOrigRcut = 1; |
| 522 |
– |
origRcut = theRcut; |
| 523 |
– |
} |
| 524 |
– |
|
| 533 |
|
rCut = theRcut; |
| 534 |
|
checkCutOffs(); |
| 535 |
|
} |
| 536 |
|
|
| 537 |
< |
void SimInfo::setEcr( double theEcr ){ |
| 537 |
> |
void SimInfo::setDefaultRcut( double theRcut ){ |
| 538 |
|
|
| 539 |
< |
if( !haveOrigEcr ){ |
| 540 |
< |
haveOrigEcr = 1; |
| 541 |
< |
origEcr = theEcr; |
| 534 |
< |
} |
| 539 |
> |
haveOrigRcut = 1; |
| 540 |
> |
origRcut = theRcut; |
| 541 |
> |
rCut = theRcut; |
| 542 |
|
|
| 543 |
+ |
( rCut > ecr )? rList = rCut + 1.0: rList = ecr + 1.0; |
| 544 |
+ |
|
| 545 |
+ |
notifyFortranCutOffs( &rCut, &rList, &ecr, &est ); |
| 546 |
+ |
} |
| 547 |
+ |
|
| 548 |
+ |
void SimInfo::setEcr( double theEcr ){ |
| 549 |
+ |
|
| 550 |
|
ecr = theEcr; |
| 551 |
|
checkCutOffs(); |
| 552 |
|
} |
| 553 |
|
|
| 554 |
+ |
void SimInfo::setDefaultEcr( double theEcr ){ |
| 555 |
+ |
|
| 556 |
+ |
haveOrigEcr = 1; |
| 557 |
+ |
origEcr = theEcr; |
| 558 |
+ |
|
| 559 |
+ |
( rCut > ecr )? rList = rCut + 1.0: rList = ecr + 1.0; |
| 560 |
+ |
|
| 561 |
+ |
ecr = theEcr; |
| 562 |
+ |
|
| 563 |
+ |
notifyFortranCutOffs( &rCut, &rList, &ecr, &est ); |
| 564 |
+ |
} |
| 565 |
+ |
|
| 566 |
|
void SimInfo::setEcr( double theEcr, double theEst ){ |
| 567 |
|
|
| 568 |
|
est = theEst; |
| 569 |
|
setEcr( theEcr ); |
| 570 |
|
} |
| 571 |
|
|
| 572 |
+ |
void SimInfo::setDefaultEcr( double theEcr, double theEst ){ |
| 573 |
|
|
| 574 |
+ |
est = theEst; |
| 575 |
+ |
setDefaultEcr( theEcr ); |
| 576 |
+ |
} |
| 577 |
+ |
|
| 578 |
+ |
|
| 579 |
|
void SimInfo::checkCutOffs( void ){ |
| 580 |
|
|
| 581 |
|
int cutChanged = 0; |
| 642 |
|
( rCut > ecr )? rList = rCut + 1.0: rList = ecr + 1.0; |
| 643 |
|
|
| 644 |
|
if( cutChanged ){ |
| 613 |
– |
|
| 645 |
|
notifyFortranCutOffs( &rCut, &rList, &ecr, &est ); |
| 646 |
|
} |
| 647 |
|
|