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 777 by mmeineke, Wed Sep 17 14:22:15 2003 UTC vs.
Revision 778 by mmeineke, Fri Sep 19 20:00:27 2003 UTC

# Line 290 | Line 290 | template<typename T> void Integrator<T>::moveA(void){
290    int i, j;
291    DirectionalAtom* dAtom;
292    double Tb[3], ji[3];
293 –  double A[3][3], I[3][3];
294 –  double angle;
293    double vel[3], pos[3], frc[3];
294    double mass;
295  
# Line 326 | Line 324 | template<typename T> void Integrator<T>::moveA(void){
324  
325        for (j = 0; j < 3; j++)
326          ji[j] += (dt2 * Tb[j]) * eConvert;
329 –
330 –      // use the angular velocities to propagate the rotation matrix a
331 –      // full time step
332 –
333 –      dAtom->getA(A);
334 –      dAtom->getI(I);
335 –
336 –      // rotate about the x-axis      
337 –      angle = dt2 * ji[0] / I[0][0];
338 –      this->rotate(1, 2, angle, ji, A);
339 –
340 –      // rotate about the y-axis
341 –      angle = dt2 * ji[1] / I[1][1];
342 –      this->rotate(2, 0, angle, ji, A);
327  
328 <      // rotate about the z-axis
345 <      angle = dt * ji[2] / I[2][2];
346 <      this->rotate(0, 1, angle, ji, A);
347 <
348 <      // rotate about the y-axis
349 <      angle = dt2 * ji[1] / I[1][1];
350 <      this->rotate(2, 0, angle, ji, A);
351 <
352 <      // rotate about the x-axis
353 <      angle = dt2 * ji[0] / I[0][0];
354 <      this->rotate(1, 2, angle, ji, A);
328 >      this->rotationPropagation( dAtom, ji );
329  
330        dAtom->setJ(ji);
357 –      dAtom->setA(A);
331      }
332    }
333  
# Line 666 | Line 639 | template<typename T> void Integrator<T>::constrainB(vo
639      painCave.isFatal = 1;
640      simError();
641    }
642 + }
643 +
644 + template<typename T> void Integrator<T>::rotationPropagation
645 + ( DirectionalAtom* dAtom, double ji[3] ){
646 +
647 +  double angle;
648 +  double A[3][3], I[3][3];
649 +
650 +  // use the angular velocities to propagate the rotation matrix a
651 +  // full time step
652 +
653 +  dAtom->getA(A);
654 +  dAtom->getI(I);
655 +  
656 +  // rotate about the x-axis      
657 +  angle = dt2 * ji[0] / I[0][0];
658 +  this->rotate( 1, 2, angle, ji, A );
659 +  
660 +  // rotate about the y-axis
661 +  angle = dt2 * ji[1] / I[1][1];
662 +  this->rotate( 2, 0, angle, ji, A );
663 +  
664 +  // rotate about the z-axis
665 +  angle = dt * ji[2] / I[2][2];
666 +  this->rotate( 0, 1, angle, ji, A);
667 +  
668 +  // rotate about the y-axis
669 +  angle = dt2 * ji[1] / I[1][1];
670 +  this->rotate( 2, 0, angle, ji, A );
671 +  
672 +  // rotate about the x-axis
673 +  angle = dt2 * ji[0] / I[0][0];
674 +  this->rotate( 1, 2, angle, ji, A );
675 +  
676 +  dAtom->setA( A  );    
677   }
678  
679   template<typename T> void Integrator<T>::rotate(int axes1, int axes2,

Diff Legend

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