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

Comparing trunk/OOPSE/libmdtools/Integrator.cpp (file contents):
Revision 1097 by gezelter, Mon Apr 12 20:32:20 2004 UTC vs.
Revision 1125 by gezelter, Mon Apr 19 22:13:01 2004 UTC

# Line 182 | Line 182 | template<typename T> void Integrator<T>::integrate(voi
182    // initialize the forces before the first step
183  
184    calcForce(1, 1);
185
186  //temp test
187  tStats->getPotential();
185    
186    if (nConstrained){
187      preMove();
# Line 213 | Line 210 | template<typename T> void Integrator<T>::integrate(voi
210    MPIcheckPoint();
211   #endif // is_mpi
212  
213 <  while (info->getTime() < runTime){
213 >  while (info->getTime() < runTime && !stopIntegrator()){
214      if ((info->getTime() + dt) >= currStatus){
215        calcPot = 1;
216        calcStress = 1;
# Line 690 | Line 687 | template<typename T> void Integrator<T>::rotationPropa
687  
688    double angle;
689    double A[3][3], I[3][3];
690 +  int i, j, k;
691  
692    // use the angular velocities to propagate the rotation matrix a
693    // full time step
# Line 697 | Line 695 | template<typename T> void Integrator<T>::rotationPropa
695    sd->getA(A);
696    sd->getI(I);
697  
698 <  // rotate about the x-axis
699 <  angle = dt2 * ji[0] / I[0][0];
700 <  this->rotate( 1, 2, angle, ji, A );
701 <
702 <  // rotate about the y-axis
703 <  angle = dt2 * ji[1] / I[1][1];
704 <  this->rotate( 2, 0, angle, ji, A );
698 >  if (sd->isLinear()) {
699 >    i = sd->linearAxis();
700 >    j = (i+1)%3;
701 >    k = (i+2)%3;
702 >    
703 >    angle = dt2 * ji[j] / I[j][j];
704 >    this->rotate( k, i, angle, ji, A );
705  
706 <  // rotate about the z-axis
707 <  angle = dt * ji[2] / I[2][2];
710 <  this->rotate( 0, 1, angle, ji, A);
706 >    angle = dt * ji[k] / I[k][k];
707 >    this->rotate( i, j, angle, ji, A);
708  
709 <  // rotate about the y-axis
710 <  angle = dt2 * ji[1] / I[1][1];
714 <  this->rotate( 2, 0, angle, ji, A );
709 >    angle = dt2 * ji[j] / I[j][j];
710 >    this->rotate( k, i, angle, ji, A );
711  
712 <  // rotate about the x-axis
713 <  angle = dt2 * ji[0] / I[0][0];
714 <  this->rotate( 1, 2, angle, ji, A );
715 <
712 >  } else {
713 >    // rotate about the x-axis
714 >    angle = dt2 * ji[0] / I[0][0];
715 >    this->rotate( 1, 2, angle, ji, A );
716 >    
717 >    // rotate about the y-axis
718 >    angle = dt2 * ji[1] / I[1][1];
719 >    this->rotate( 2, 0, angle, ji, A );
720 >    
721 >    // rotate about the z-axis
722 >    angle = dt * ji[2] / I[2][2];
723 >    this->rotate( 0, 1, angle, ji, A);
724 >    
725 >    // rotate about the y-axis
726 >    angle = dt2 * ji[1] / I[1][1];
727 >    this->rotate( 2, 0, angle, ji, A );
728 >    
729 >    // rotate about the x-axis
730 >    angle = dt2 * ji[0] / I[0][0];
731 >    this->rotate( 1, 2, angle, ji, A );
732 >    
733 >  }
734    sd->setA( A  );
735   }
736  

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines